{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 11.1基本的采样方法"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "在预测分布中，后验概率组要用于计算期望，我们希望找到某个函数$f(z)$关于一个概率分布$p(z)$的期望\n",
    "$$\\mathbb E[f]=\\int f(z)p(z)dz$$\n",
    "简单来说采样方法就是得到概率分布$p(z)$中独立抽取的一组变量$z$,其中$l=1,...L$,这样期望就可以通过有限和的方式得到:\n",
    "$$\\hat f=\\frac{1}{L}\\sum_{l=1}^Lf(z^{(l)})$$\n",
    "只要样本$z^{(l)}$是从$p(z)$中抽取的，那么$\\mathbb E[\\hat f]=\\mathbb E[f]$,因此估计$\\hat f$有正确的均值，方差就为\n",
    "$$vat[\\hat f]=\\frac{1}{L}\\mathbb E[f-\\mathbb E[f])^2]$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "最早的蒙特卡罗方法是为了求解一些积分问题"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "import sys\n",
    "#sys.path\n",
    "sys.path.append(r\"../\")\n",
    "import matplotlib.animation as animation\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "%matplotlib inline\n",
    "from prml.rv import Gaussian, Uniform\n",
    "from prml.sampling import metropolis, metropolis_hastings, rejection_sampling, sir\n",
    "\n",
    "np.random.seed(666)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.text.Text at 0x7f38584e4390>"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAPUAAAD+CAYAAADidkOmAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAE2hJREFUeJzt3XtsVOWfx/HPM21DJK6Ayj9uMCBZBSNBV8XNxhjcxGi8\n/NXN+lPZzWJiXP1HNzFeolbcX9w/TH76j3FXJD9dKSg/0R8Gud9aLpVelBbojV6hpffbtJ3pZeY8\n3/3jmVJqC53bmTmXzyshpsNMe0r79jlzznOeo0QEROQdgWxvABGlF6Mm8hhGTeQxjJrIYxg1kccw\naiKPYdREHsOoiTyGURN5TG4iT37iiSdk//79dm0LEV2biveJCUXd19eX+KYQ0SwiAkyMAyPDkNFh\nYHQEMjoChEOQ0CgwFoKEQ0A4BIyNAaGRv8v74E8N8XzuhKImovhJOAQZ6AMG+iBDA5ChASA4CBka\nAkaGgEhk7heq2KA8dV2G+bgOQE48X5dRE6VIwiFIdwekuxPS0wXp7Qb6uoHxsZlPvDrWQADIieWn\nremAp/5+xhcQAAjHuz2MmigBMjkB6WiHtLdCOtogl9uA4aHpJygFCMw74JxcAAJYVuzFV8WqNQBt\nyzYyaqLrkEgEcqkZ0twAfbEJ6GyPBQkz2mpt4hWZOeIKACualW1m1ES/I8FB6PpqSEMtpKUhNtIq\nM/oqZSK2otNxZynea2HURIiFfO4MdG0V0NFuHgwEzMg7FbEgNhLbs9ucLoyafEsmJyA1Z6EryyAX\nm82DgQAQyAFEO3Ykng+jJt+R/l7ospPQleXA5EQs5NiorJ09CseDUZNv6EvN0CePQhpqp08v5eRM\nH532CEZNnqdbGqGL9kMutZgRWQXM7jXguaABRk0epjvaoA/9DGltNO+Tp05BeRyjJs+R4SCsI3sg\nZ3+9KmbvjcjXwqjJM0Rb0KUnoYv2m3nVSvkq5imMmjxBujsR/elboPPy9PtmHwYNMGpyOdEa+tRR\n6KIDM2d7+RijJteS4SFYP24zE0emDoLxNlKMmtxJN9bB+nGbWWjAg+eaU8GoyVVEBPrkEeij+0zM\nSjHo32HU5BoSicD66VtIdZXZ3WbMc2LU5AoSGoX13Z8h7Rd9M4kkWYyaHE+GBhDd+gUwOMD3z3Fg\n1ORo0ttlgg6FgADfP8eDUZNjSVcHot/8j7k8Uinfn3+OF+/QQY4k3VcFDTDoBDBqchzp7Ub0m/+9\nKmjucieCUZOjSHAQ0cIvzKQSgEEngVGTY0h41BwUGx0BwINiyWLU5AgSjcD69s/AQD8PiqWIUVPW\niWhYu76LTSzhCJ0qRk1Zp4sOQqorOfUzTRg1ZZWuPQt9/FBsUQNO/UwHRk1ZI71dsHZ9a6Z+CoNO\nF0ZNWSET44ju+NrsbnNdg7Ri1JRxIgJrz06gv9c84NO1xOzCqCnj5EwZ5NwZHhizCaOmjJLeLlj7\n/sprom3EqCljxIoi+uN2xmwzRk0Zo4sOAl2Xzc3bGbZtGDVlhG5rgT511JyP5vtoWzFqsp1EIrB2\nfRe7BzRHaLsxarKdLj4ADPTB7HeT3Rg12Uo62qBLinjlVQYxarKNaAvR3X8xu92UMfzXJtvo0pNA\nVwcAxXtcZRCjJlvI8JC5T3QgwN3uDGPUZAvrwE/mxu88OJZxjJrSTjdfgNScNQfHeLFGxjFqSivR\nFqz9u8w10pw1lhWMmtJKl5cAvd3gbnf2MGpKGwmHeHDMARg1pY0uPgSMj5v53ZQ1/NentJCBPuiK\nU7ELNjhKZxOjprSwjuzJ9iZQDKOmlOm2VnMKC+BVWA7AqCklIgJ9eA9PYTkIo6aUSFM95FIz53Y7\nCKOmpIloWEf2AgGO0k7CqClpUntues0xcgxGTUkRrWEd229Gaa455iiMmpIi1ZVAX4+5aIMchVFT\nwkRrWMUHY6M0J5o4DaOmhMn5M+Y+WAGO0k7EqCkhojWs44fMKB3lKO1EjJoSIjVVHKUdjlFT3EQ0\nrBNHOEo7HKOmuMmFGqCnk+elHY5RU1xEBNbxwzwv7QKMmuIiLQ1AR1u2N4PiwKgpLvrU0dgcb47S\nTseoaV7S0QZpbgDAK7HcgFHTvKySIrOYIK/EcoXcbG8AOZsM9Jlz0+QaHKnT6Ouvv8aaNWuwcOFC\nKKVQV1eHoqIi5OXloaWlJeHPt337dixatAj9/f02bG189C/FWfvalBxGnSZ79+7Fxo0bsWbNGnz+\n+efYunUr7rzzTrz11lt44YUXsGLFioQ/57PPPoulS5fio48+smGL5yfhEHRleewDvp92C+5+p8mX\nX36J22+/HYWFhQjE7sd85MgRlJWV4bPPPkvqc+bk5ODll1/GBx98gIKCAixevDidmzwvXVECRCNm\n/TGem3YNjtRpUlJSgocffvhK0IAJfeXKlXjwwQeT/rx/+MMfMD4+jm3btqVjM+Mm0Qh02cnYOt4M\n2k0YdYpee+01KKXQ09OD7du3QymF3NxcjIyMYPfu3Xj88cdnPH9wcBCLFy/G+vXrZzze3d2NFStW\nYPXq1RgYGLjy+LJly7B69Wr88MMPmfh2rpBzvwGhUXPUm1yFP7EU5efn48033wQAvP3229i6dSt2\n7NiB6upqhMNhrFu3bsbzlyxZgjfeeAPFxcUoKioCAIRCITz11FMYHx/Hvn37cPPNN894zUMPPYTS\n0lJEIpGMfE8iAuv0cS6C4FKMOkWPPPIIbrvtNgDAq6++ig0bNiA/Px+1tbUAgDvuuGPWa15//XUs\nXboUmzZtgmVZePbZZ1FfX4+9e/di+fLls56/cuVKhMNhtLa22vmtXCEtDUBPFy/ccClGnQZVVVVY\nsmQJli1bduWx3t5eAGZk/r0bb7wR77zzDoqLi/HYY4/hwIED2LlzJ+677745P/8tt9wCAOjr67Nh\n62fTp4/H7lzJ99JuxKjToLKyEmvXrk3oNa+88goWLVqEY8eOYfPmzbPee19NZ3Aml/T3QhpqM/b1\nKP0YdYoikQhqampmRb106VIA5sDYXD755BMEg0EAwE033XTdrzH1OW699dZUN3deuvSEWSGUU0Jd\ni1GnqK6uDhMTE7j33ntnPL569WoAQFNT06zXFBYW4t1330VBQQHuuusuFBQUXHc0bmpqwsKFC+d8\nv51OMj42PdmEXItRp6iyshIAZo3U999/P2644QaUl8+M5PDhw3jxxRexceNGfPjhh3j//fdRU1OD\nwsLCa36N0tJSrFu3Dnl5een/Bq6iz5QBkUlz1Jtci1GnqLKyErm5ubj77rtnPJ6Xl4dnnnkGBw8e\nvPJYVVUV8vPzsX79emzevBkA8Nxzz2HVqlXYtGnTnKes2traUFtbi/z8fFu/D9Eauvxk7AAZT2O5\nGaNOUVVVFVatWoUFCxbM+ruXXnoJjY2NqKioQFtbG5588kksX74cO3fuRG6umaEbCARQUFCAlpYW\nbNmyZdbn2LFjBxYsWIANGzbY+n1IQy0wOGBmkJGrKUlgov4DDzwgFRUVNm6Ot4gI1q1bh3vuuQdf\nffVVwq+3LAurVq3C008/jU8//dSGLZwW3foFpLWJK5s412jeB3/6m3ieyP8t20gphY8//hiFhYVJ\nXXq5Y8cO9PT04L333rNh66ZJbzek+QK4sok3MGqbPfroo4hEIkldevn8888jGAxemXxiF11+KnYa\ni1F7AaP2OZkYh66aOkLPqL2AUfucrqoAJid5NZaH8CfpYyJidr05z9tTGLWPSWsTbxzvQYzaxzhK\nexPXKPMpGR6C1J3L9maQDThS+5T+9bRZIZS73p7DqH1ILMtErXjXDS9i1D4kdeeB0AiQwx+/F/Gn\n6kPmAFkOEOXVWF7EqH1GersgF5vA2WPexah9RpeXcLkij2PUPiKTE2ZaKHkao/YROXcGmJzgckUe\nx6h9QkRgVZzickU+wKh9Qi5fAro6ONnEBxi1T+iKEt7B0icYtQ9IOAQ5fwY8jeUPjNoHdGWZGaG5\nUqgv8KfscSIauuIXc4CMK4X6AqP2OGluAAb7OUr7CH/SHqcrSngay2e4SIKHyfAQpL4625tBGcaR\n2sOuLIQAnpv2E0btUWJZ0L/FFkIQXrzhJ4zao6T+PDDKhRD8iD9xj9LlJVwIwacYtQdJbzektRGc\nQeZPjNqDpm94x/fSfsSoPYYLIRCj9hh99rfYQgj80foVf/IeIiLQFbyVjt8xag+RSy1AdycXQvA5\nRu0huvwkR2li1F4hI0FI7bnYtFDyM0btEfrX0+YUFne9fY9Re4BYllkIgTe8IzBqT5Das7zhHV3B\n3wIP0KUnOM+brmDULqc72iDtF7O9GeQgjNrldNnJ2HtpnsYig1G7mIRGuJ43zcKoXUxXnOZ63jQL\nfxtcSqLR2CWW3PWmmRi1S0l1JU9j0Zz4G+FCIgLr9HGexqI5MWoXkkstQNdlrvxLc2LULqRPH+fV\nWHRNjNplZKAPUsersejaGLXL6NPHzZVYbJqugVG7iIyFzb2mzUdZ3RZyLkbtIrqiBIhEuKggXRd/\nO1xCohHo0tg8bx4go+tg1C4hVb+aySYcpWke/A1xAdEaVskxM9mEN4+neTBqF5C688BAH9cfo7gw\naocTEVinjnKUprgxaoeTlgagoy3bm0EuwqgdTp84YkZpXl5JcWLUDqbbWnmfaUoYo3YwfeKwOYXF\ntbwpAYzaoaTrMqShNtubQS7EqB3KKj7EUZqSwqgdSLoum8sriZLAqB3IKj7IUZqSxqgdRjrbzQwy\noiQxaofhKE2pYtQOottaIfXV2d4McjlG7RAiAn10b2z2GEdpSh6jdghpvgBpbQJnj1GqGLUDiGhY\nR/cCORylKXWM2gGkugroaM/2ZpBHMOosk2gU1pE9seuleSUWpY5RZ5kuPwUMDXJVE0obRp1FMhaG\nPh6b481VTShNGHUW6WP7gfEx3jSe0oq/TVkiPZ1mcX7FUZrSi1FngYjA2r8rdk8snpem9GLUWSB1\n5yAtjVMfZXVbyHsYdYbJ5IQZpTnRhGzCqDNMFx0AhoMAeAqL7MGoM0i6O6bvL82DY2QTRp0hojWs\nn3eao918G002YtQZoktPQNovxva6WTXZh1FngAz0QR/dF5s5xvndZC9GbTMRDWv397xtDmUMo7aZ\nLj0Zu3WO4iksyghGbSPp6YI+vCe2kCBHasoMRm0TsaKI/nWbmQbK42KUQYzaJvrwHqCrwxztFu52\nU+Ywahvo+vNmkgmPdlMWMOo0k+AgrF3fcW43ZQ2jTiOJRmD95f+AyYlsbwr5GKNOExGBtecHSEdb\nbG43d7spOxh1mujyU5DK8thKJgyasodRp4FurIPev8scGOORbsoyRp0i6emE9f03vFMlOQajToEM\nDyG6fYu5NpprjZFDMOokSTiE6NYvgJFhc2CMozQ5BKNOgkyMw9q+BejvNUFHuYoJOQejTpBMTsDa\nvsWcuuKMMXIgRp0AmZyAtW0LpK2VQZNj5WZ7A9xCxsJmhG6/aKaAMmhyKEYdBxkdRnTrZqC3i0GT\n4zHqeUhvlzltNRzkLje5AqO+Dt18wVygEY1wPje5BqOeg4hAnz4OfWh3bOonuBwRuQaj/h2ZnID1\n8/eQc2diu9saXI+I3IRRX0U62xH9odBMKuFcbnIpRg1zSxz9SzH0sX3mAR7hJhfzfdTS1w3rpx3m\n/HMgEFv9k7vb5F6+jVoiEegTh6FLjpm3zBydySN8F7WIQGrOwjq8GxgaNCuVAHz/TJ7hq6h1axP0\n0b1m7nZODpCTy/tEk+f4Imp9sRn6+EFIcwMQyOGuNnmaZ6MWrSEXaqBLjk1fVcV7WpEPeC5qCY1C\nV5VDl58y75kDOeYPYyaf8ETUYkUhjXXQVb9C6s+bg16BAJCby1VJyHdcG7VYUUhzA3TtWUjtOWB8\nbPo889TIzCPa5EOuilqCg5DmC9ANdZCmenN7G6XMX149Kgt3tcm/HBu1iADBQUhbK+RiM3Rro5mT\nDUyPyFeHzN1sIgAOiVpEgNFhSFcHpLPd/Gm7CIRGzBOmRuOcHEDL9G41QyaaJaNRi7aA4BCkvxfS\n3wv09UB6uyE9ncBYePqJU1dI5eaa88lTc7F5bploXmmLWrQGxkLAyAhkJAiMBCHBIcjwkNmNHhwA\ngoMzL5aYGoEBE7DW038AjsRESUgs6tAIrKP7gHAIMhYGxkKQ0REgNGpG2rmubpoKd+qodI4yCw+I\nnvl8BkyUFglFLcEg9InDM0MFzMeBgPmvYOZ86qvD1RbAs0xEtkpu9/v3I7II3+8S2Ssn3icqSWBB\ngIV5edbKWxZHktokIkraWCQijf1DN8Tz3ISiVkpViMgDSW8ZESUlkfZ4Ly0ij2HURB6TaNSbbdkK\nIppP3O0l9J6aiJyPu99EHsOoiTyGURN5TFxRK6X+WSl1WilVpZQ6qZRaaveGERGglFqvlNqayGvi\nHamPicg/iMhaAIcA/EvCW0dEyVgL4EwiL4g36n9XSpUppaoAvApgPNEtI6Kk3Avgb5VSpUqpZqXU\n+vleMG/USql/A7AOwD/FRup6ANWpbikRxWUtgBEReQjAfwD443wviGekXgOgRERGlVL5AP4RwLmU\nNpOI5qWUygNwK4D/jj1UGfv4uuKJ+msAryqlygDcB6BZREJJbicRxW8VgEYRmYx9/PcAquZ7EWeU\nETmUUupfAfwXTNx5MAep/1NETl/vdY5YTZSI5rQWwI8ASgDcAOCP8wUNcKQm8hzOKCPyGEZN5DGM\nmshjGDWRxzBqIo9h1EQew6iJPIZRE3nM/wPZ/IxlbjKX4AAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f385b86b4a8>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "def sigmoid(x):\n",
    "    return np.exp(x) / (1 + np.exp(x))\n",
    "x=np.linspace(-5,5,100)\n",
    "y=sigmoid(x)\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(4, 4))\n",
    "ax.spines['right'].set_color('none')\n",
    "ax.spines['top'].set_color('none')\n",
    "ax.set_xticks([-5,5])\n",
    "ax.set_xticklabels(['$a$','$b$'])\n",
    "ax.set_yticks([])\n",
    "ax.plot(x,y,color=\"salmon\")\n",
    "ax.fill_between(x,y,color=\"salmon\")\n",
    "ax.set_xlim(-5,5)\n",
    "ax.set_ylim(0,1)\n",
    "ax.text(-2, 0.6, \"$f(x)$\", fontsize=\"xx-large\")\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "则一个简单的近似求解方法是在[a,b]之间随机的采样一个点。比如 $x_{0}$ ，然后用 $f(x_{0})$ 代表在[a,b]区间上所有的 $f(x)$ 的值。那么上面的定积分的近似求解为:\n",
    "$$(b-a)f(x_{0})$$\n",
    "推广开来，我们可以采样[a,b]区间的n个值： $x_{0}，x_{1},...,x_{n-1}$ ，用它们的均值来代表$[a,b]$区间上所有的 $f(x)$ 的值。这样我们上面的定积分的近似求解为:\n",
    "$$\\frac{b-a}{n}\\sum_{n=0}^{n-1}f(x_i)$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "但是绝大部分情况下x在区间$[a,b]$的分布不是均匀分布，如果假设x的分布为$p(x)$则：\n",
    "$$\\int_a^bf(x)dx=\\int_a^b\\frac{f(x)}{p(x)}p(x)\\simeq\\frac{1}{n}\\sum_{n=0}^{n-1}\\frac{f(x_i)}{p(x_i)}$$\n",
    "代入均匀分布页可以得到前面的结果"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "现在考虑如何从简单的分均匀分布产生随机数，我们使用一个函数$y=f(z)$来对z服从的$（0,1）$均匀分布变换\n",
    "$$p(y)=p(z)|\\frac{dz}{dy}|$$\n",
    "$$z=h(y)\\equiv\\int_{-\\infty}^yp(\\hat y)d\\hat y$$\n",
    "常用的函数包括\n",
    "$$p(y)=\\lambda exp(-\\lambda y)$$\n",
    "$$p(y)=\\frac{1}{\\pi}\\frac{1}{1+y^2}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "例子:使用Box-muller方法用于生成高斯分布的样本"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYEAAAD8CAYAAACRkhiPAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAFXlJREFUeJzt3X+MndWd3/H3Jy4laBMUEF7X2KYmklPJuBunjFxX9A8K\nzcbNRjGpVOSoG1gV4UiwCVEjJTiRmqSVJaTdJLuohcrZIMwuG6+lJMKiUNawrKJIa5wx4ZcNNFaB\n4pHB3s1GBlVy1863f9zDcpmMmTs/78w875d0Neee+zzPPQOe+7nPOec5T6oKSVI3vWfYDZAkDY8h\nIEkdZghIUocZApLUYYaAJHWYISBJHWYISFKHGQKS1GGGgCR12D8YdgMmc8kll9TatWuH3QxJWlQO\nHTr011W1fLLtFnwIrF27ltHR0WE3Q5IWlSSvDLKd3UGS1GGGgCR1mCEgSR1mCEhShxkCktRhhoAk\ndZghIEkdNnAIJFmW5KdJHmzPL06yP8nP2s+L+rbdkeRokheTfKyv/sokz7bX7kyS2f11JElTMZUz\ngduA5/ue3w48VlXrgMfac5KsB7YBVwBbgLuSLGv73A3cDKxrjy0zar0kaUYGumI4yWrgt4CdwH9s\n1VuBq1t5N/CXwJdb/Z6qOg28lOQosCnJy8CFVXWgHfM+4Drg4dn4RaS5dsdP//odz2//yCVDaok0\newZdNuIPgC8B7++rW1FVx1v5NWBFK68CDvRtd6zV/V0rj6+XFiVDQUvBpN1BST4BnKiqQ+fapqoK\nqNlqVJLtSUaTjJ48eXK2DitJGmeQMYGrgE+27pw9wDVJ/gR4PclKgPbzRNt+DFjTt//qVjfWyuPr\nf0VV7aqqkaoaWb580kXwJEnTNGkIVNWOqlpdVWvpDfj+RVX9NrAPuLFtdiPwQCvvA7YlOT/J5fQG\ngA+2rqNTSTa3WUE39O0jSRqCmSwlfQewN8lNwCvA9QBVdTjJXuAIcAa4tarOtn1uAe4FLqA3IOyg\nsJYMxwi0GKXXnb9wjYyMlPcT0DCM/1CfKUNB8ynJoaoamWw7rxiWpA4zBCSpwwwBSeowQ0CSOswQ\nkKQOm8kUUWnJme0ZQdJC55mAJHWYISBJHWYISFKHGQKS1GEODEvzpH/Q2SUktFAYAuo0ZwOp6+wO\nkqQOMwQkqcMMAUnqMMcEpCHwBjRaKAwBdYoDwdI7TdodlOS9SQ4meTrJ4STfaPVfTzKW5Kn2+Hjf\nPjuSHE3yYpKP9dVfmeTZ9tqd7V7DkqQhGeRM4DRwTVW9meQ84MdJ3ro38Ler6vf7N06ynt4N6a8A\nLgUeTfKhdp/hu4GbgSeAh4AteJ9hSRqaSc8EqufN9vS89ni3GxNvBfZU1emqegk4CmxKshK4sKoO\nVO/GxvcB182s+ZKkmRhodlCSZUmeAk4A+6vqifbS55I8k+SeJBe1ulXAq327H2t1q1p5fP1E77c9\nyWiS0ZMnT07h15EkTcVAIVBVZ6tqI7Ca3rf6DfS6dj4IbASOA9+crUZV1a6qGqmqkeXLl8/WYSVJ\n40xpdlBV/SLJ48CW/rGAJN8BHmxPx4A1fbutbnVjrTy+XpozzgaS3t0gs4OWJ/lAK18AfBR4ofXx\nv+VTwHOtvA/YluT8JJcD64CDVXUcOJVkc5sVdAPwwCz+LpKkKRrkTGAlsDvJMnqhsbeqHkzyx0k2\n0hskfhn4LEBVHU6yFzgCnAFubTODAG4B7gUuoDcryJlBkjRE6U3UWbhGRkZqdHR02M3QIrVYu4O8\nglgzleRQVY1Mtp1rB0lShxkCktRhhoAkdZghIEkdZghIUocZApLUYYaAJHWYISBJHWYISFKHeXtJ\naQHyHsSaL4aAlpTFukyENCx2B0lShxkCktRhhoAkdZghIEkdZghIUocNcnvJ9yY5mOTpJIeTfKPV\nX5xkf5KftZ8X9e2zI8nRJC8m+Vhf/ZVJnm2v3dluMylJGpJBzgROA9dU1YeBjcCWJJuB24HHqmod\n8Fh7TpL1wDbgCmALcFe7NSXA3cDN9O47vK69LkkakkmvE6je/SffbE/Pa48CtgJXt/rdwF8CX271\ne6rqNPBSkqPApiQvAxdW1QGAJPcB1+F9hqVJefGY5spAYwJJliV5CjgB7K+qJ4AVVXW8bfIasKKV\nVwGv9u1+rNWtauXx9ZKkIRkoBKrqbFVtBFbT+1a/YdzrRe/sYFYk2Z5kNMnoyZMnZ+uwkqRxprRs\nRFX9Isnj9PryX0+ysqqOJ1lJ7ywBYAxY07fb6lY31srj6yd6n13ALoCRkZFZCxctTS4VIU3fILOD\nlif5QCtfAHwUeAHYB9zYNrsReKCV9wHbkpyf5HJ6A8AHW9fRqSSb26ygG/r2kSQNwSBnAiuB3W2G\nz3uAvVX1YJK/AvYmuQl4BbgeoKoOJ9kLHAHOALdW1dl2rFuAe4EL6A0IOygsSUM0yOygZ4CPTFD/\nN8C159hnJ7BzgvpRYMOv7iFJGgavGJakDjMEJKnDDAFJ6jBDQJI6zNtLSotQ/7URLiGhmfBMQJI6\nzBCQpA4zBCSpwxwT0KLjWkHS7PFMQJI6zBCQpA4zBCSpwwwBSeowQ0CSOszZQVrwnA0kzR3PBCSp\nwwwBSeqwQe4xvCbJ40mOJDmc5LZW//UkY0meao+P9+2zI8nRJC8m+Vhf/ZVJnm2v3dnuNSxJGpJB\nxgTOAF+sqieTvB84lGR/e+3bVfX7/RsnWQ9sA64ALgUeTfKhdp/hu4GbgSeAh4AteJ9haUbGj5m4\nqqimYtIzgao6XlVPtvIbwPPAqnfZZSuwp6pOV9VLwFFgU5KVwIVVdaCqCrgPuG7Gv4EkadqmNCaQ\nZC29m84/0ao+l+SZJPckuajVrQJe7dvtWKtb1crj6yd6n+1JRpOMnjx5cipNlCRNwcAhkOR9wPeB\nL1TVKXpdOx8ENgLHgW/OVqOqaldVjVTVyPLly2frsJKkcQYKgSTn0QuA+6vqBwBV9XpVna2qXwLf\nATa1zceANX27r251Y608vl6SNCSDzA4K8F3g+ar6Vl/9yr7NPgU818r7gG1Jzk9yObAOOFhVx4FT\nSTa3Y94APDBLv4ckaRoGmR10FfAZ4NkkT7W6rwCfTrIRKOBl4LMAVXU4yV7gCL2ZRbe2mUEAtwD3\nAhfQmxXkzCBJGqJJQ6CqfgxMNJ//oXfZZyewc4L6UWDDVBooSZo7XjEsSR3mAnJacFwwTpo/hoC0\nxHgFsabC7iBJ6jBDQJI6zBCQpA4zBCSpwwwBSeowQ0CSOswQkKQOMwQkqcMMAUnqMENAkjrMEJCk\nDjMEJKnDDAFJ6rBJVxFNsga4D1hB7y5iu6rqD5NcDPwZsJbencWur6q/bfvsAG4CzgKfr6pHWv2V\nvH1nsYeA26qqZvdXktTPVUX1bgZZSvoM8MWqejLJ+4FDSfYDvwM8VlV3JLkduB34cpL1wDbgCuBS\n4NEkH2q3mLwbuBl4gl4IbMFbTHae9w+QhmfS7qCqOl5VT7byG8DzwCpgK7C7bbYbuK6VtwJ7qup0\nVb0EHAU2tRvTX1hVB9q3//v69pEkDcGUxgSSrAU+Qu+b/IqqOt5eeo1edxH0AuLVvt2OtbpVrTy+\nXpI0JAOHQJL3Ad8HvlBVp/pfa9/sZ61vP8n2JKNJRk+ePDlbh5UkjTNQCCQ5j14A3F9VP2jVr7cu\nHtrPE61+DFjTt/vqVjfWyuPrf0VV7aqqkaoaWb58+aC/iyRpiiYNgSQBvgs8X1Xf6ntpH3BjK98I\nPNBXvy3J+UkuB9YBB1vX0akkm9sxb+jbR5I0BIPMDroK+AzwbJKnWt1XgDuAvUluAl4BrgeoqsNJ\n9gJH6M0surXNDAK4hbeniD6MM4MkaagmDYGq+jGQc7x87Tn22QnsnKB+FNgwlQZKkuaOVwxLUocZ\nApLUYYaAJHWYISBJHTbI7CBJS4gLyqmfZwKS1GGGgCR1mN1BGgqXj5YWBs8EJKnDDAFJ6jBDQJI6\nzBCQpA4zBCSpwwwBSeowQ0CSOszrBKSO679mwyUkusczAUnqsEnPBJLcA3wCOFFVG1rd14GbgZNt\ns69U1UPttR3ATcBZ4PNV9Uirv5K3by35EHBbVdVs/jJauLxCWFqYBjkTuBfYMkH9t6tqY3u8FQDr\ngW3AFW2fu5Isa9vfTS841rXHRMeUJM2jSUOgqn4E/HzA420F9lTV6ap6CTgKbEqyEriwqg60b//3\nAddNt9GSpNkxkzGBzyV5Jsk9SS5qdauAV/u2OdbqVrXy+HpJ0hBNNwTuBj4IbASOA9+ctRYBSbYn\nGU0yevLkycl3kCRNy7RCoKper6qzVfVL4DvApvbSGLCmb9PVrW6slcfXn+v4u6pqpKpGli9fPp0m\nSpIGMK0QaH38b/kU8Fwr7wO2JTk/yeX0BoAPVtVx4FSSzUkC3AA8MIN2S5JmwSBTRL8HXA1ckuQY\n8DXg6iQbgQJeBj4LUFWHk+wFjgBngFur6mw71C28PUX04faQJA3RpCFQVZ+eoPq777L9TmDnBPWj\nwIYptU6LltcFSIuDy0ZI+nvjw9tlJJY+l42QpA4zBCSpwwwBSeowQ0CSOsyBYUnn5EDx0ueZgCR1\nmCEgSR1mCEhShzkmoFnhFcLS4uSZgCR1mCEgSR1mCEhShxkCktRhhoAkdZghIEkdZghIUodNGgJJ\n7klyIslzfXUXJ9mf5Gft50V9r+1IcjTJi0k+1ld/ZZJn22t3tnsNS5KGaJAzgXuBLePqbgceq6p1\nwGPtOUnWA9uAK9o+dyVZ1va5G7iZ3s3n101wTEnSPJs0BKrqR8DPx1VvBXa38m7gur76PVV1uqpe\nAo4Cm5KsBC6sqgNVVcB9fftIkoZkustGrKiq4638GrCilVcBB/q2O9bq/q6Vx9dPKMl2YDvAZZdd\nNs0mSpptLi299Mx47aCqqiQ1G43pO+YuYBfAyMjIrB5bs8O1gqSlYbqzg15vXTy0nyda/Riwpm+7\n1a1urJXH10uShmi6IbAPuLGVbwQe6KvfluT8JJfTGwA+2LqOTiXZ3GYF3dC3jyRpSCbtDkryPeBq\n4JIkx4CvAXcAe5PcBLwCXA9QVYeT7AWOAGeAW6vqbDvULfRmGl0APNwekqQhmjQEqurT53jp2nNs\nvxPYOUH9KLBhSq2TJM0pbyojadqcLbT4uWyEJHWYISBJHWZ3kAbmtQHS0uOZgCR1mCEgSR1mCEhS\nhxkCktRhDgzrnBwIlpY+zwQkqcM8E5A0a/rPHr16eHHwTECSOswQkKQOMwQkqcMMAUnqMAeGJc0J\nl5leHGYUAkleBt4AzgJnqmokycXAnwFrgZeB66vqb9v2O4Cb2vafr6pHZvL+ml1eFyB1z2x0B/2r\nqtpYVSPt+e3AY1W1DnisPSfJemAbcAWwBbgrybJZeH9J0jTNxZjAVmB3K+8Gruur31NVp6vqJeAo\nsGkO3l+SNKCZhkABjyY5lGR7q1tRVcdb+TVgRSuvAl7t2/dYq5MkDclMB4b/ZVWNJfl1YH+SF/pf\nrKpKUlM9aAuU7QCXXXbZDJsoSTqXGYVAVY21nyeS/JBe987rSVZW1fEkK4ETbfMxYE3f7qtb3UTH\n3QXsAhgZGZlyiGgwDgRLmnYIJPk14D1V9UYr/ybwn4F9wI3AHe3nA22XfcCfJvkWcCmwDjg4g7ZL\nWkScMrowzeRMYAXwwyRvHedPq+p/JvkJsDfJTcArwPUAVXU4yV7gCHAGuLWqzs6o9ZKkGZl2CFTV\n/wY+PEH93wDXnmOfncDO6b6nJGl2ecVwhzgGoIXE7qGFwbWDJKnDDAFJ6jBDQJI6zDGBJc5xAC0W\njhEMh2cCktRhhoAkdZghIEkd5pjAEuMYgJaK/n/Ljg/MHUNA0oLnoPHcsTtIkjrMM4FFzu4fSTPh\nmYAkdZhnAouM3/wlxwhmkyGwwPmhL03OUJg+u4MkqcMMAUnqsHnvDkqyBfhDYBnwR1V1x3y3YaGx\ny0eaXe/2N2VX0TvNawgkWQb8N+CjwDHgJ0n2VdWR+WzHfPNDXlo4HD94p/k+E9gEHG33JybJHmAr\nvZvPLxl+6EuLx2R/r0s9JOY7BFYBr/Y9Pwb883luw0D8IJcEU/ssWIyBsSCniCbZDmxvT99M8uI0\nD3UJsBA/zW3X1NiuqbFdUzNr7doxGwd520zb9Y8H2Wi+Q2AMWNP3fHWre4eq2gXsmumbJRmtqpGZ\nHme22a6psV1TY7umpuvtmu8poj8B1iW5PMk/BLYB++a5DZKkZl7PBKrqTJLfBR6hN0X0nqo6PJ9t\nkCS9bd7HBKrqIeCheXq7GXcpzRHbNTW2a2ps19R0ul2pqvl4H0nSAuSyEZLUYZ0JgSRfTFJJFsRE\n3iT/JckzSZ5K8udJLh12mwCS/F6SF1rbfpjkA8NuE0CSf5fkcJJfJhn6TI4kW5K8mORoktuH3R6A\nJPckOZHkuWG3pV+SNUkeT3Kk/T+8bdhtAkjy3iQHkzzd2vWNYbepX5JlSX6a5MG5fJ9OhECSNcBv\nAv9n2G3p83tV9RtVtRF4EPhPw25Qsx/YUFW/AfwvZn3q87Q9B/xb4EfDbkjf8if/BlgPfDrJ+uG2\nCoB7gS3DbsQEzgBfrKr1wGbg1gXy3+s0cE1VfRjYCGxJsnnIbep3G/D8XL9JJ0IA+DbwJWDBDIBU\n1am+p7/GAmlbVf15VZ1pTw/Qu5Zj6Krq+aqa7kWDs+3vlz+pqv8HvLX8yVBV1Y+Anw+7HeNV1fGq\nerKV36D3wbZquK2C6nmzPT2vPRbE32GS1cBvAX801++15EMgyVZgrKqeHnZbxkuyM8mrwL9n4ZwJ\n9PsPwMPDbsQCNNHyJ0P/UFsMkqwFPgI8MdyW9LQul6eAE8D+qloQ7QL+gN4X11/O9RstyGUjpirJ\no8A/muClrwJfodcVNO/erV1V9UBVfRX4apIdwO8CX1sI7WrbfJXeafz989GmQdulxSvJ+4DvA18Y\ndyY8NFV1FtjYxr5+mGRDVQ11TCXJJ4ATVXUoydVz/X5LIgSq6l9PVJ/knwKXA08ngV7XxpNJNlXV\na8Nq1wTup3ftxLyEwGTtSvI7wCeAa2se5xBP4b/XsA20/IneluQ8egFwf1X9YNjtGa+qfpHkcXpj\nKsMeWL8K+GSSjwPvBS5M8idV9dtz8WZLujuoqp6tql+vqrVVtZbeafs/m48AmEySdX1PtwIvDKst\n/dpNf74EfLKq/u+w27NAufzJFKT3Dey7wPNV9a1ht+ctSZa/NfstyQX07nMy9L/DqtpRVavbZ9Y2\n4C/mKgBgiYfAAndHkueSPEOvu2pBTJsD/ivwfmB/m77634fdIIAkn0pyDPgXwP9I8siw2tIGzt9a\n/uR5YO9CWP4kyfeAvwL+SZJjSW4adpuaq4DPANe0f1NPtW+5w7YSeLz9Df6E3pjAnE7HXIi8YliS\nOswzAUnqMENAkjrMEJCkDjMEJKnDDAFJ6jBDQJI6zBCQpA4zBCSpw/4/+arDx3xcrIcAAAAASUVO\nRK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f3858207c18>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "def getNormal(SampleSize):\n",
    "    iid = np.random.uniform(0,1,SampleSize)\n",
    "    normal1 = np.cos(2*np.pi*iid[0:int(SampleSize/2-1)])*np.sqrt(-2*np.log(iid[int(SampleSize/2):SampleSize-1]))\n",
    "    normal2 = np.sin(2*np.pi*iid[0:int(SampleSize/2-1)])*np.sqrt(-2*np.log(iid[int(SampleSize/2):SampleSize-1]))\n",
    "    return np.hstack((normal1,normal2))\n",
    "SampleSize = 100000\n",
    "normal = getNormal(SampleSize)\n",
    "plt.hist(normal,np.linspace(-4,4,81),histtype='bar',facecolor=\"skyblue\")\n",
    "plt.show()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 11.1.2拒绝采样"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "设定一个程序可采样的分布 $q(x)$ 比如高斯分布，然后按照一定的方法拒绝某些样本，以达到接近 $p(x)$ 分布的目的，其中$q(x)$叫做 proposal distribution。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "具体算法如下：\n",
    "* 引入一个常数k\n",
    "* 从概率分布$q(z)$中产生一个数$z_0$\n",
    "* 在区间$[0,kq(z_0)]$上的均匀分布生成一个数$\\mu_0$\n",
    "* 如果$\\mu_0>\\tilde p(z_0)$则样本被拒绝"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAW4AAAD8CAYAAABXe05zAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xl4W+Wd9//3fSR53/c1cRInzr5AoCSBkAKllLZAm1KW\nQqFlCuWhT0s7ZWY6nT5dZrrMNZ0W2umPGQYYmJY2kFIKBRogJJBAyGaI4+zxnjje93iVzrl/f8g2\ncWLHsi3pSM73dV26bMtH53xsS18f3edelNYaIYQQ4cOwO4AQQoiJkcIthBBhRgq3EEKEGSncQggR\nZqRwCyFEmJHCLYQQYUYKtxBChBkp3EIIEWakcAshRJhxBmKnaWlpuqCgIBC7FkKIaam4uLhZa53u\ny7YBKdwFBQXs3bs3ELsWQohpSSlV7eu20lQihBBhRgq3EEKEGSncQggRZqRwCyFEmJHCLYQQYUYK\ntxBChJmAdAcUQoSHzs5OGhsbcbvddkeZ1lwuFxkZGSQkJPhlf1K4RXir2Qm170P2MshZARExdicK\nG52dnTQ0NJCbm0t0dDRKKbsjTUtaa3p7e6mtrQXwS/GWwi3CU20xbPkxlL85fJdWBu6cS3Dd8hQq\nIcfGcOGhsbGR3NxcYmLkn10gKaWIiYkhNzeXU6dOSeEWF6g3vg/vPkyfEcs7xjpKrblk0sRM1cDK\n2r30/MdaYu7fjEousDtpSHO73URHR9sd44IRHR3ttyYpKdwivHzwO3j3YfYZS9ik19GPCxR0Es9x\nPZtDFHDHwAv0/uZKou59HSOjyO7EIU2aR4LHn79r6VUiwseJPeiXv0m1UcBL1tXeon2WUyqbp7gZ\n7Rmg/7+uxtNea0NQIQJLCrcID5116GfvoIs4NvJJtBr7qduo0nmaz+Eye2j47b1BDClEcEjhFuHh\nxQcwe9p4Vt1Et44cd/MmlcYOVpLb8g5tH/wlCAGFCB4p3CL0VW6H8jd5S63mlJnk88O2cyltJMIr\n30J7+gMYUISCu+66i+TkZB555BH+8Ic/cPfdd4/4vtaa5cuX8/TTT/u8z6997Wvcc889fk46dT4X\nbqWUQyn1gVLq5UAGEmIErWHLP9PtSGSXuXhCD/UoF3/loyR7Gqn70/cCFFCEgra2Nl577TUef/xx\nfve73/HAAw/w5S9/ecQ2zz33HK2trdx+++0+7/fb3/42zzzzDGVlZf6OPCUTOeP+BnA4UEGEGNXx\n1+HELt7mMjzq3IuR4z5czeYIc0g79CS9jZUBCChCQXJyMvX19axfv549e/bQ2trK2rVrR2zzq1/9\nijvvvBOXy/fnUUFBAZdffjmPPvqovyNPiU+FWymVB3wSeDywcYQ4g2Wh3/wRHY5Uis35k97N61yJ\nCzd1L/7Qj+FEKElJSeF73/seX/jCF0hNTSUlJYX77ruP/n5vE1lZWRk7duzgc5/73IjHPfvssyil\nRr11dnYCsH79ep555hksywr6zzUWX8+4Hwb+Dgid5GL6O/QCquEAW1mNpRyT3k2bSuIIc8iq3YS7\np8OPAUUoqKmpoa2tjV/+8pekpKSwceNGHnzwQR577DF+85vfAPDmm28SGxvLsmXLRjx27dq1vPfe\ne8O3V155hZSUFK6//vrhEY6rV6+moaGB0tLSoP9sYxl3AI5S6lNAo9a6WCm17jzb3QvcCzBjxgy/\nBRQXKK3R235OqyODErMQpjh2YRcXsYByqv76MAXrv++fjCIklJSUAN4LiT/72c8AuOqqq9i2bRub\nNm3iW9/6FsXFxSxYsADDGHmump2dTXZ2NuCdu+Xqq69m8eLFbNy4cXibRYsW4XA42L179zmF3y6+\njJxcA9yglLoeiAISlFK/01rfceZGWuvHgMcAVq5cqf2eVFxYTu5BNR5il/M68MOIs2ryqCOd+IO/\nRX/meyhDOlSN5od/OcihU522HHthTgLf//SiCT+upKSEmJgYHnrooRH3z5s3j927dwNQX19PWlra\nmPvo7u7m+uuvxzAMXn755RHztzidTpKSkqivr59wtkAZ99mrtf6O1jpPa10A3ApsObtoC+F3xU/j\nVpGUeAr9sz+l2MVFpFpNnHr39/7ZpwgJJSUlrF27ltTU1BH319fXk5PjnWysr6+PyMjR+//39fVx\nww030N3dzaZNm4iPjz9nm8jISPr6+vwffpJkrhIRevo60Af+SClFDKgIv+32AEVcw3asd/8DrpBz\nj9FM5ozXbiUlJVxxxRUj7uvq6mLLli388IfeC9IpKSmjnjG73W7Wr19PXV0db7/9NsnJyaMeo729\nnZSUFP+Hn6QJvV/UWr+ltf5UoMIIAUDpRpSnj32Gf9sTTeVkL8vI7ztM6/Hdft23sEd3dzfl5eWU\nl5ePuP9Xv/oV4B2UA1BUVERl5cjuoKZpctttt3H06FE2b95Menr6qMdoamqip6eHefPmBeAnmBxp\n6BOhRWt08VM0GlmcMMduk5ysYpZgoWje8v/5fd8i+EpLS7Esi/Lycn70ox+xdetWvvOd7/DDH/6Q\nxx57jKQk70jbNWvWUFNTQ1NT0/Bj77//fl599VV+8pOfUFNTw86dO9m5c+c5vUf27t2LUorVq1cH\n9Wc7H2kqEaHl1Aeo+lLed14L2v9Tjp5WcVTpfNLrt2KZJoZj8t0Mhf1KSkpISEjglVde4Utf+hI/\n/vGPKSoqYuPGjdx4443D261bt46UlBQ2bdrEnXfeidaaDRs20Nvbyy233DJin5/+9Kd56aWXhr/e\ntGkTV1555Tlt6HaSwi1CS/FTeFQE+zxzp9wFcCylzOdG/Ton975M3kduHP8BImSVlJSwZMkSli5d\nSnFx8ZjbRUREcMcdd7BhwwbuvPPOEQNszsc0TZ5//vnhboahQppKROjw9KMPvsAh5tKvxp8BcLIO\nU4gHB727fxuwY4jgKCkp8blv9UMPPcTWrVs5duyYz/vfuHEj0dHR3HrrrZONGBBSuEXoqHgb1d/J\nEcfCgB6mX0VxnFlkt+zAMyCzBoYrrTWlpaUsXbrUp+3z8vJ48sknqaurm9AxnnjiCZzO0GqcCK00\n4sJ26EUGjGiOebID1kwypJT5LKCM6h3PMnPdFwN7MBEQvjZ3nGmiZ8633XbbhLYPFjnjFqHBdKOP\nvMwRPQtTBf584jiz6CcCd7EMxhHhRwq3CA2V21B97Rx1LAjK4TzKxWEKyet6n77TMvGUCC9SuEVo\nOPQibhXJUU9u0A55gPlE0c+p7XKRUoQXKdzCfqYHfeQVjhGcZpIhleTTRySeAy+Nv7EQIUQKt7Bf\nzQ5UTzNHHZNfLGEyLOXgOAXkdJfilt4lIoxI4Rb2O/QiHhXBYU9+0A99jNnE0cOpPbKUqggfUriF\nvSwLffgvHKdgUmtKTlUZs7BQ9Jb8KejHFmKypHALe9XvR51uoNxpz8xrfSqKGnJJadqN1rL+hwgP\nUriFvcreAOCIO3i9Sc52jNlk6EYaju6xLYMQEyGFW9hKH99MvZFNt4q1LcMxZgPQsec52zIIMRFS\nuIV9elrh5G7KmGVrjBaVQgtJRNVstTWHEL6Swi3sU7EVpS3KlL2FG+AYc8h1V9LRVGt3FCHGJYVb\n2Of4ZvqMGGrM0ZeMCqajzMaJSdNOaS4JRxs2bCAyMhK32z2l/WitWb58OU8//bRP23/ta1/jnnvu\nmdIxJ0MKt7CHZaGPv0GZnolW9j8Na8iljwjMo6/ZHUVMQklJCQsXLsTlmlqX0ueee47W1lZuv/12\nn7b/9re/zTPPPENZWdmUjjtR9r9ixIWpvgTV00Slc47dSQDQyqCKfDJOH8I0TbvjiAkqKSlhxYoV\nU97Pr371K+68806f/wEUFBRw+eWX8+ijj0752BMhhVvY4/hmNIoj7jy7kwyrYCbJdFB/6D27o4gJ\n2rdvH8uXLx/+esuWLaSmpvLggw9imiYdHR3cd999pKSkkJSUxLe+9S1++tOfUlhYOPyYsrIyduzY\nwec+97kR+3722WdRSo166+zsZP369TzzzDNYlhW0n1cWUhC20Mdfp97IpkfH2B1lWDkzAejc9yK5\nSy63OY3wVVNTE3V1dcOF+9e//jV///d/zyOPPMJXvvIV3G431113Hd3d3fznf/4nsbGx/NM//RMN\nDQ2sWrVqeD9vvvkmsbGx5yyFtnbtWt5778N/5q2trdx5551cdtllJCQksHr1ahoaGigtLfV5GbWp\nksItgq+3HWr3Uq4uszvJCK0k0U48ESfesTuKff76D1Bfas+xs5bAJya+KG9JSQlKKRYuXMhXvvIV\n/vSnP/Hqq6+ybt06AH7+859TVlbG0aNHSUlJASA5OZk1a9aMWPasuLiYBQsWYBgjGyKys7PJzs4G\noLOzk6uvvprFixezceNGABYtWoTD4WD37t1SuMU0VvUOSltUGDMhlEaZK0WFnsnCgeP0dncRHRtv\ndyLhg3379pGSksJnP/tZmpqa2LVr13ATiGVZPPzwwzzwwAPDRRtg9mzvoKszC219fT1paWljHqe7\nu5vrr78ewzB4+eWXiYnxvlt0Op0kJSVRX18fiB9vVFK4RfBVvo1bRVBtZgR8bcmJqmAmF3GAyvf/\nyqwrPm93nOCbxBmv3UpKStBas337djZv3jyi3frAgQM0NjZyzTXXjHjMyZMngZGFu6+vb7gYn62v\nr48bbriB7u5utmzZQnz8yH/qkZGR9PX1+etHGpdcnBRBpyvepoZcLOWwO8o5KpiBBvoO/tXuKMJH\nJSUlfP3rX+eOO+7gS1/6Ek1NTcPfq631DqjKyxt5EXzz5s0kJCRQUFAwfF9KSgrt7e3n7N/tdrN+\n/Xrq6up4/fXXSU5OPmeb9vb2EWf0gSaFWwRXZx2q+ShVRoHdSUbVq6KpI5O4xr12RxE+6O/v58iR\nIyxdupTHHnuMtLQ0br75ZjweDwCpqakAVFZWDj+mubmZX/ziFyxZsgSlPnzLV1RUNGI7ANM0ue22\n2zh69CibN28mPf3cwWJNTU309PQwb17wZriUwi2Cq3IbABUEf9EEX1UwgxzrJK31NXZHEeM4dOgQ\nbrebJUuWEB0dzQsvvMDBgwf55je/CXibQnJycnjwwQd55ZVXePbZZ7n22msxTfOcC4lr1qyhpqZm\nxBn7/fffz6uvvspPfvITampq2LlzJzt37qS09MMLuHv37kUpxerVq4PzQyOFWwRb5dv0GTGcMlPt\nTjKmcmbiwKKl+M92RxHjKCkpITY2ljlzvAO5Zs6cybPPPsujjz7KU089RWRkJH/84x8BuPnmm/nl\nL3/Jww8/jGmaI3qUAKxbt46UlBQ2bdoEeIe/b9iwgd7eXm655RZWrVo1fPvud787/LhNmzZx5ZVX\nDp/dB4MUbhE8WqMr3qJC54EKsauSZzhBDm6ceI69aXcUMY67776b06dPj2jyuOqqq/B4PNx9990A\nrFq1ipKSEnp6eti5cyf5+fl0dHScc8YdERHBHXfcwYYNGwCGB9horc+5vfSSd4Fp0zR5/vnngz5f\niRRuETytFajOWmoc9s8GeD6mcnKCHFI6D8qqONPQvn37MAyDJUuWnPO9hx56iK1bt3Ls2DGf9rVx\n40aio6O59dZb/R3zvKRwi+CpeAuAMtO+1W58VUUembqJxuojdkcRflZSUsLs2bOJjT138Y68vDye\nfPJJ6urqfNqX1ponnngCpzO4PaulH7cInsq36TISabESQ67/9tmqBi+etu97hcyCBTanEf70gx/8\ngB/84Adjfn8iZ8+33XabHxJNnJxxi+CwLHTFNip0fki3bw+pJYsBnFgVb9sdRYhzSOEWwdF4ENXX\nRo2zwO4kPrGUgxPkktp1KKizvgnhCyncIjiqvBM3lXuybQ7iuyryydDNNFUdtjuKECNI4RbBUfUO\nHUYyHYTPxE1VeIdJt+972eYkgSO9ZoLHn7/rcQu3UipKKbVbKVWilDqolPqh344uLgyWha7eQaXO\nsTvJhJwikwFcWIOjPacbl8tFb2+v3TEuGL29vVNeWm2IL2fc/cBVWutlwHLgOqVCbCJlEdqaDqN6\nWzkZJu3bQyzloJpc0roOT8t27oyMDGpra+np6ZEz7wDSWtPT00NtbS0ZGRl+2ee43QG19y96evBL\n1+BN/srCd1XvAuHVvj2kmjzm8g71FQfIKlw6/gPCSEJCAgCnTp2a8uro4vxcLheZmZnDv/Op8qkf\nt1LKARQDhcBvtNa7/HJ0cWGo2k6nkUS79s+TNtBOmfG8785lgbORHIe3P3dHycvTrnCDt3j7q5iI\n4PGpcGutTWC5UioJeEEptVhrfeDMbZRS9wL3AsyYMcPvQUWY0jps2rd7tZPd7nwqzFQUmh3umWQa\nXfSrCHTldrvjCTFsQr1KtNbtwFbgulG+95jWeqXWeuVoc9aKC1TTUVRPM7WOAruTnFe/dvDnvkVU\nmclcFNnAI9em4nQYvO0upIYcUk4fm5bt3CI8+dKrJH3wTBulVDTwMUAmcBC+qfKeqZaZod2+XW6m\n0oeLL+a18b8PfZ4brlrFjz+zhEYrll3mfDJoprFKnvYiNPhyxp0NbFVK7Qf2AG9oradvx1bhX9Xv\n0mUk0hbC7dtaw3FPGmmOHu658aPExcUB8NmLZ3D13CT+7LkUgI5SWc5MhIZxC7fWer/WeoXWeqnW\nerHW+kfBCCamAa3Rle9QpXNCen6SFh1Dq45hSUzXiLUJlVL8+22XUO6cQ792YZZPz/7cIvzIyEkR\nOC1lqJ4mToZ4+/ZxTxoOLG5bM3fEhPwASTERXL2kgPetQpI6D0l/ZxESpHCLwKke7L9tZtkcZGwe\nragwUyhwtrN65YpRt/nYoix26QVk6SZaaitH3UaIYJLCLQKnegc9RhwtOsnuJGOqNpMZwMnaXOdw\n2/bZrpibRrFegKE0rfulnVvYTwq3CBhdvSPk27ePm2nEq35uueqiMbeJiXASMeMS3NrBwLG3ghdO\niDFI4RaB0X4C1XGCWsdMu5OM6bTlos5KoCiynXlz5553248um81+PZvYNmnnFvaTwi0Co+Y9ACqt\n0O2/XW95p5i9dmEmhnH+l8LVCzLYbc0njzraGk8FI54QY5LCLQKj+l36VRT1VordScbUYMUTgcma\nxbPH3TY7MZra+GW4lElLyWtBSCfE2KRwi4DQ1TuoIQetQvcp1mDFkek4zawC35pzspdehakVzSVy\ngVLYK3RfVSJ8nW5CNR/jpBG6k431aScdOpr86AFiYmJ8eszaJYUc1jNJOF0W4HRCnJ8UbuF/g+3b\n1eTaHGRsDZa369/FM3zvqrg4N4H9xgLmcILOtuZARRNiXFK4hf/VvIdHuThhptmdZEwNVhwGFlcu\nmeXzY5RSdGZcQpRyU/+BtHML+0jhFn6nq9/lJNlYymF3lDE1mnGkO3oonDWx7orRc68AoL74lUDE\nEsInUriFf/V1Qn0pJ1S+3UnG5NYGzTqGHFcviYmJE3rsvNlzKLeySTx9PEDphBifFG7hXyd2o7TF\nCSNv/G1t0mzFojFYmhVzzqRS41mSl8geaz6zOUFvT3eAEgpxflK4hX9Vv4uFQZXHP6tZB4L3wqTm\no0smPqozLtJJZcxi4lUv9aVv+T2bEL6Qwi38Ste8R72RhVu57I4ypgYrjhSjjwWFBZN6vJ65BoDO\n/Zv8mEoI30nhFv7j7oPa4pDuBmhpaLTiyHKcJiNjcu8K8mcv4JROwTq5x8/phPCNFG7hP7XFKHMg\npAfetOoYPDiYnxox4fbtIcvzk9ljzaeAkwz09/s5oRDjk8It/KdmBwAVIdy+3Wp5R0leNm/yk18V\nZcXzvi4iWXXRcGSXv6IJ4TMp3MJ/qnfQZGTQp6LtTjKmFisGFybLCiffnBPhNGhOuRiAroOv+yua\nED6Twi38w/Sga3ZRrXPsTnJerTqaFEcvOdlTm242Y85y2nQcVuU7fkomhO+kcAv/qN+PcndT6wzd\nhRO09jaVpDn6xlymzFfLZqSwxyoic6AK0zT9lFAI30jhFv4xOLFUuTt027e7dCQeHBQkOSd9YXLI\n0rxEdlnzSVdtNJSV+CmhEL6Rwi38o3oHHUYyXSre7iRjatXeC5PLZqROeV8FqbEccCwEoLNU+nOL\n4JLCLaZOa3T1e1SGePt2ixWNQnNx4dRzGobCkbOMbh2JrnrXD+mE8J0UbjF1zcdQvS2cCuGFgcHb\nvp1s9DEj1z/rYM7PTaXYmkfy6SNYluWXfQrhCyncYuqqvWecFWaWzUHOr9WKIVn1kJo69aYSgPnZ\n8eyyFpBFM001x/yyTyF8IYVbTF3Vu3Qb8bToiU2RGkx92kkPEeTFMe6K7r6anxXPbms+AO0yb4kI\nIincYmq0Rle9Q6XOhSn21AikVss7KGhhdoLf9jk3I55SZjOgnZgV2/y2XyHGI4VbTE1bJep0PbUh\n3r7dMjjUfWWhf9q3AaIjHGQkxXNIF5DUcQittd/2LcT5SOEWUzPYo6Lc9F9BDIRWHUOsGmDeTP/2\nfFmQk8Quaz5Zup7Wumq/7luIsUjhFlNTvYNeI5YmnWx3kvNqtWJIMXrJzMz0634XZiey3VqMgaal\nRNq5RXBI4RZToqveoUrnhHT7tkcrOnQUWZFuIiIi/Lrvoqx43rfmYWqFu+wtv+5biLFI4RaT134C\n1VET8u3b7ToajaIwzf+zFi7IjqeHKCrJJaHtgLRzi6CQwi0mr3pw/m0rtEdMDvUouWiW/+dRyU+O\nIcqpKKGIbOsU7U11fj+GEGeTwi0mr/od+lU09VaK3UnOq92KxonFktn+v4BqGIq5GXG8ay7AiSnt\n3CIopHCLSdNVO6gmB61C+2nUpqNJNvrI8vOFySELc5N4y1yChWLg2JsBOYYQZwrtV5wIXV31qNYy\nToTw+pJD2qxokoxeEhMDM7JzQVY8rTqeWrKIb5EpXkXgjVu4lVL5SqmtSqlDSqmDSqlvBCOYCHGD\n85OE8oruAH3aQS8R5MQy5Tm4x1KU5R2NeUjP8rZzN9cH5DhCDPHljNsD/K3WeiFwGfCAUmphYGOJ\nkFe5nQEVyUkzze4k59U+eGFyXmbg5gmfn+Xd9/vMx4lJ8wevBuxYQoAPhVtrXae1fn/w8y7gMIT4\naZYIOF21nRpyw6J9G2BFAHqUDEmOjSAt1skuayEWiv5jWwJ2LCFggm3cSqkCYAWwa5Tv3auU2quU\n2tvU1OSfdCI0ddahWsqoNkK7/zZ4z7gjlEnRjMBOOTs/K4EaM4k6MkiQdm4RYD4XbqVUHPA88KDW\nuvPs72utH9Nar9Rar0xPT/dnRhFqqrYDUEm+zUHG16ajSVb+H+p+tgU5iXRYUVSRR7ZVK+3cIqB8\nKtxKKRfeov2M1vpPgY0kQl7lNvpVNKfM0O6/rbW3R0mKo4/Y2NiAHmtuZjwmBof0bGnnFgHnS68S\nBTwBHNZa/yLwkUSo01XbqQqD9u1eXAzgJC/eEbAeJUOGLn7u03OlnVsEnC+vvDXAncBVSql9g7fr\nA5xLhKr2E6i2Kk44CuxOMq62wR4l8/24eMJY5mbEAdCgUqWdWwScc7wNtNbvAKE79ZsIrsH27XIr\n9DsWDRXuFbMDvxZmbKST7IQI2nqjqXbkcam1j/amOpLSQ3uechGeQvu9rgg9ldvpM2JCfn4S8M4K\nGKPcFOYHZxHj+VkJtHoiqWAmTkyail8KynHFhUcKt/Cd1ujKbVTqvJCef3tImxVNkuolWL2cirIT\n6dDeniUmBp6jbwTluOLCI4Vb+K6tCtV5Mizat7WGdh1FmmuA6Gj/z8M9mnmZcVgoWnQcJ8gmub1U\n5ucWASGFW/husH27zAzt+bcBTusIPDjITxz3Mo7fDPUsadPRVDKDTF1P84ljQTu+uHBI4Ra+K9/K\naSMh5NeXhA+Hui/KDV7Wwow4DAWdKo4KZqKA1uIXg3Z8ceGQwi18Y1noirep0Plh0r4dA8DyWYEd\nMXmmKJeDnMRI2sxITpFJPy50xVtBO764cEjhFr6p34/qbaHaOdvuJD5p09HEGwMU5AW3O96C7ERa\nzUgs5aCaPNK7DmFZVlAziOlPCrfwTcVWAI66Q799G7w9SpJVT9B6lAwpykqgU0dhakUlM0iljfqj\ne4OaQUx/UriFb8q30mRk0K0CO+eHP5ha0aEjSY/w4HK5gnrsuZlxaBQdOooKvKsDde6T/tzCv6Rw\ni/G5e9E1OynXoT8bIECHjkJjMCslMujHLhpcVKHdiqaRNE4Tg6P6naDnENObFG4xvpr3UGY/VY5Z\ndifxydBQ96Uzgr86z6y0WG/PEiMWlKKSfLL6jjHQ3x/0LGL6ksItxle+FVM5qfCEx7wbbToaA4tl\nQZij5GyRTgczk6NoM6MAKKeAeLqpK9kc9Cxi+pLCLcalK7ZykhzcKrjtxZPlXdW9n9zs4BdugPk5\nibSY3maacryrBPWUSDu38B8p3OL8Tjeh6kupNArsTuKzNu2doyQlxZ6JsOZnJdKlI3Frg9MqjnrS\niavfaUsWMT1J4RbnV/k2AOWDPSRC3YB20K0jyYnRGIY9T2/vBUpF++DozTIKyDFraGs4YUseMf1I\n4RbnV7aZPiOGWjP4F/omo83yti0XpsfYlmFB9mDPEu3NUEYBDixa9siqf8I/pHCLsVkW+vgblOmZ\nIb9M2ZChOUqWz8qwLUN+cgxRTkU73j7vJ8ihHxfm0ddsyySml/B4NQp71H2A6mmmwjnX7iQ+a7Ni\ncGGyqMC+HjCGoShMjx2+QGkph3e2wK4DmB6PbbnE9CGFW4zt2OtoFEfcob9M2ZA2HU2y0UdmZvAm\nlxrNorxkWq1ohqbjLqeAJDqoOyiDccTUSeEWY9LHX+eUkUOvCs5CBFOltXfEYoqjl/j4eFuzzM+K\np1876cXbhbJssFvg6Q9esDOWmCakcIvRnW5EnXqf4yo8ZgME6MVFP07y4gyUzVPPzs/yriw/NIqz\nXSXRTDJRJ+WMW0ydFG4xujLvSL/jhMcwd/iwSM7PTrQ5ifeMG6DjjEm5yikg11NFV2ujXbHENCGF\nW4zu+Ot0G/GcMlPtTuKz1sHCvWK2ve3bAMmxEaTGOIcXdAA4xmxceGjc+ayNycR0IIVbnMv0oMve\n5JguCIu4tgySAAAauUlEQVTVboa06hhi1QBzZ4TGnCoLshOGe5YAVJFHHxGYh162MZWYDqRwi3Od\n2IXq76TCUWh3kglpsWJINXrIyrJnjpKzLcxNot2KwhrsWWIpB2XMIuf0ftz9ffaGE2FNCrc41/HX\nsHBwzBMeq90AuLVBh44iO8pDRESE3XEAKMqMx8SgU0cN33eUOcTRQ+2ev9iYTIQ7KdxiJK3Rh1+m\nWuUzoIK/EMFkedu3FUUZ9g11P9v8waHvQxdNAY5TgIlB3z4Z/i4mTwq3GKnpCKq1nKOOIruTTEjL\n4LwglxSGRvs2QGFGHIaCDiNu+L5+FUU1eaS17JJFhMWkSeEWIx3+CxrFIbPA7iQT0mrFEqU8LJoV\nOs07kU4HM1OiaDWjRtx/lNmk6RYaDu2wKZkId1K4xQj68F+oNXLpIvQXBT6T98JkN9nZoXPGDbAw\nJ4mWcwr3HAA69mywI5KYBqRwiw+1VaHq93NUhc+kUuBd1b1NR5EZ4SY6OrSG5y/KTeS0jqRfO4bv\n61CJ1JNG3MltNiYT4UwKt/jQYW//4kN6js1BJqZNR6MxmJsaehdTl+R6R3G2WCMvmh5lDjlmDa0n\nj9kRS4Q5KdziQ0deptHIpFUn2J1kQoaK4so5odF/+0yLc7yFu42Rk14dZh4GmpbtT9mQSoQ7KdzC\nq6sBXbOTI2HWTALewh2BybI5oXNhckhybATZCRE065Fn3A2keSedKn/VpmQinEnhFl5HX0GhOUx4\njZaEwQuTjh5yckKvcAMsy0+myXNW27tSHKSIXE8VrTVH7AkmwpYUbuF16EXajVTqzWS7k0yIpaFN\nx5Du7CcuLm78B9hgaX4SXWddoAQ4ONRcsu0Jm5KJcCWFW0BXPbpyG6XMC6tJpQA6dBQmBnOSXXZH\nGdNYFyibVBqNpBJdKWtRiomRwi3gwPMobVGqFtqdZMJaLG9/84tmpducZGxDhbuVc1flOUgRuWY1\nzRX7gx1LhLFxC7dS6kmlVKNS6kAwAong0/ufo97IpsmyfwGCiWq2YnBicdHcPLujjCkpJoKchIjh\nYflnOkARCmh958ngBxNhy5cz7qeA6wKcQ9il+Tiqbh+laoHdSSalwYojw9HNjPzQLdwAS0e7QAm0\nqmTqSCeu+g300MrCQoxj3MKttd4GtAYhi7DD/ufQKEqtc7sBurVBmSeFD9zZNFsxhFpdGdAGbTqG\n3Ig+2xcHHs9YFyjB21ySY56k+fheG5KJcOT0146UUvcC9wLMmDHDX7sVgaQ1uvQ5qo2ZdOkP5yZp\nsaLZ787mhJWEOfi/fZ8nl3jVxxxHC8uc9RjK/ireZMWhUSzLDc3eJGc68wJljqNrxPdKmc/VvEP7\n24+SPu8SO+KJMOO3i5Na68e01iu11ivT00P3QpE4Q20xqq2KA8aHFyUbrVj+2l9EnU6gyNXC/XN7\neO7WAr6yLJrUCIt9nlzeGpiNqe3vfdJgxaHQrFs80+4o4zrfBcpOlUA5M8mo3YzpcQc7mghDfjvj\nFmFo/3OYyskBzyxQ0GDG8sbAPKKVh68v0dx18xdwubzd7C5dvojvAo/8tYRfvg1bB+awLqIcp41n\n3g1WHKlGL/MLQ38l+uELlL2jL/Swj0V8jlc58e4G8q+8M8jpRLiR7oAXKncvev+zHGUO/SqSRjOW\n1wfmEWN4+OZyxT23fXa4aJ/pG59Yxnc+VsAJK4ktA4W2nXmbWtFkxZLl7CY1NTxWoh/rAiXAEQrp\nJRL3nqeCG0qEJV+6A/4BeA8oUkqdVErdE/hYIuAOvYTqa+cDxwoGtIOtA3OIMTx8a7nB3Z+/CXWe\ngTj3Xb2I/3fdbGqtRPZ57Jn/ulXHYOJgYXrkebOGkuUzkunSkfTqc9/omspJKfPJP72P0821NqQT\n4cSXXiW3aa2ztdYurXWe1lrG504Hxf9Du5FKmSeb3e48enFxc95pvnjzjT4Vwi+vW8AnF6ZS6smm\nyQr+ogsNpveC5LpF+UE/9mRdUpACeJt4RvMBi3HhoWHzb4IZS4QhaSq5EDUehpr3KFZLqLUSOW6m\nsyyika9/wbeiPeSnn7+Y5CiD7e5ZeILcZNJgxZFg9LN0XkFQjzsVS3ITiXQomhh9oFM9GdSTRuzx\nF6RPtzgvKdwXouKnMJWTXeYi3nEXkGT08XefXEJi4sRGTiZEuXjk9pV0WFEUe4I3AEZrb+HONE6H\n3FJl5xPhNFien0idZ4x3KEqxj8Vkmado3P9mcMOJsCKF+0Lj7kWX/IEjFPK2ey692sVnc7tZdenK\nSe3uinkZfH5FJoc8mTSYwWky6dBR9ONibpLC4Th3QEsoW1WYTosVPepAHIASFuDGSc9bDwc5mQgn\nUrgvNAdfQPV1sE1dyjEznYWuZr7+hRumdIHvB59ZTnKUYo9nRlBGVzYOthGvmRd6K96M59JZKYAa\n/hnO1qeiKWEh+W076G6qCW44ETakcF9ItIY9T9BmpPH7vo/gwuLe1XkkJSVNabcxEU7+8VOLabJi\nqTBT/BR2bHVWPNHKzUcWhn7/7bOtyE/GaTBmOzfALlbgxKThlZ8FMZkIJ1K4LyQ170HtXt7QH+GE\nlcyKqCauu+oKv+x6/UX5zEuLpNiTh1sH7mllacVJM5F8R2dYTq0QHeFgUXb82O3cQLNKpYyZpFe9\niKe/J4jpRLiQwn0heedh+oxYHun7BLFqgP973RIiI/2zMrphKP5l/Qq6dQQHzcA1YdRZ8Qzg5KJM\nJxEREQE7TiCtLkyn2Yo57z+4nVxEPKepfUO6BopzSeG+UDQcguOv8VfrMmp1CqviWlh1ycV+PcSl\ns1K5Zl4y+91ZdOvArEhTYybhxOSmy+YHZP/BcOnsVCzUefu/l1PgXUz4g/9BW1YQ04lwIIX7QrHj\n13hUBL/ov5Fko5cHb1odkB4Z379pOSjFBwHoHqi1t3DnObpYsrDI7/sPlotnJmNw/nZulGIXF5Fp\n1tKw589ByybCgxTuC0HHSXTpc2zRKzmh01iX0sWihYFZOCE/JYYvXjaD454UWq3R5+WYrGYdSw8R\nLE42iY0N/mhNf0mIcjEvI4Y68/zT0e5jIaeJwdr60yAlE+FCCveFYOejaMvi3/s/Q5rRwwM3XRHQ\n+T0e/NgCYl2KYo9/Lx7WmEkoNDdeMsev+7XD6rkZNJox552ky6NcvMsl5PQdo2nvi0FMJ0KdFO7p\nrqsBvfcJdrKUozqXj6Z1M2dOYAtfYoyLb1xTxEkznlozwW/7rTaTyHacZuXS8FvU+GyXz03DxKDO\nOv/KPXtZymli8LzxoyAlE+FACvd09/a/oj0D/MvArWQ4urn/prVBmU3vrjWzyIx1sNeTj+WHQTkd\nViQdOpqi2D6Sk5OnvkObrZ6TRpRTUcv5p6T1KBfvcCnZ/WU07n4+SOlEqJPCPZ01l6GLn+JN61IO\nWTO4OqOP2bNnB+XQkU4H//TpJbRa0Rw306a8v2rTW6yvXxbaiwL7KsrlYO3cNKrcCeOONi1mCV3E\nYr75Y0Ju4U9hCync09mWf8ajnPzzwC1kOrq578a1QT38p5blsDQ7hmJPHn2jzEHtK0vDcTONDEc3\nl1+0yI8J7fWJpTn0aBfN+vwXWj3KxXYuJbu/nMadzwYpnQhlUrinq5PFcOjPvGheQY3O4ONZfRQU\nFAQ1glKKn9+6Erd2TGn2wGormU4dxUcST5OVFX7zk4zlqqJMDAUn9PjTBLzPElpJxPXm99CegSCk\nE6FMCvd0pDVs/j69Rhw/6b+ZPEcX9332aluizMuM567L8jnmSRte/GAitIb97mySjD6++qlVYbPa\njS8SY1ysnJFItWf86XRN5eQ11pHsaaTuxR8EPpwIaVK4p6P9z0HVdn7ruZZW4vl0gSY3N9e2OA9d\nv4i0aIP3PDOxJrjgQq2VQKuO4ZK4dhYuCN/RkmO5fmku7VYUHdb4Uw8cYzbHmEVq6X/T3yIzB17I\npHBPN90t6Ne+Q62Rx78P3EiBs50v3WjP2faQmAgnP16/nDYrekJrVGoNJZ5s4tQA91+/clqdbQ/5\n2CJv088Jy4dZFZViE+twYNL8+/sDnEyEMinc081r/4jubeen7tvpx8n6uRFkZGTYnYprF2Vx/YIU\nSjw5VJu+TSPbYMXRaMVzUUwry5cuCXBCe+QmRVOUEU2N5dvvpE0ls4OV5La8Q3vJKwFOJ0KVFO7p\npOxN2L+Brazir+6lzHW2cccN19idCvBeqPzF7ZcyN8XFtoFZ4w6Hd2uDPe58opWb+65dimFM36fq\n9UvzaDBjfJ6Y6x0upY1EjJcewOptD3A6EYqm76vhQtPXgX75m7Q50vlZ/3oUmjuXJ4XUYJUol4Pf\nfXUtMS6DLQOF9I2xfJepFVsG5tCiY7gqoYmPXLwiyEmD66YVuSgUx03f3hm5lYs/8QnizHaanv5S\ngNOJUCSFezrQGl78GrSf4H+sGzhuZbIisombPxkaZ9tnykyI4n/uWUWPjuDF/kXnNJtYGt4amM0p\nK5GPxp7ix/ffjNM5+T7g4WBmaixrZidxxJN23rlLznRS5bCNj5BZv4XWbf8d4IQi1Ejhng52PgqH\nX2KzsZb/7VtFvOrnG59YQnS0f2fn85eLC1L53ZcvJt6l2DJQyJb+ORzypLNrIJ9X++dTYyWzJvoU\n//Z/1k95WbVw8ZV18+jVLqpM398hbeMyTpBNzJbvMtBUEcB0ItRI4Q53NbvQb3yPY8Y8HndfS7uO\nYV1SG6snuWp7sKyal832//dp7lqeyEkrkV3umRwz0zGVg8ujT/GLr95Aaur55/GYTq4oTCMvwcUR\nK9Pnx2hl8AKfQGHR/cRN6P7TAUwoQsn0fg863XXUojfeRadK4A/m9XzgziPH0cXf3nJ1WFzMi3A6\n+OGtl3PfNR0cOXacGRlJpKSkkJiYGJBFHkKZYSjuuXIuP/zLIZqsGNIN39aabFNJPK8/ya19L9L2\n+HpS7n8VjAvrd3chCv1Xtxhddwv6t5/B093GRnUDrw8UYWJw+zwj6EPbpyonLZGrVq+ksLCQlJSU\nC65oD/ncxXlEOeCoNbFh/cfVbDbxUVKadtK64f8EKJ0IJVK4w1F/FzyzHqulgmeNz/BmfxG1ViKr\nYxq45/OftjudmKT4KBfrL86j3J1EzwTX7NyjlrOTi0g5toGOTT8JUEIRKqRwh5uBbvjDbei6Ep53\n3MD77gJ2u/PIc3Ty/dvXhewFSeGbr6wtRCnFXvfEJ+V6nbUcZC6JO/+V06/+QKaAncakcIeT043w\n1KfQVe/yknE9BzwFbBuYhRPNV1fEUVhYaHdCMUUFabH8zZqZlJup1JnnXx3nbFoZPM8nKWEBcbt/\nSecLfyvFe5qSwh0umo6hH78Gs/4Af3TexAfmPN5xF9Cs47g6sZFbbrzO7oTCT7758YVkxnon5fK1\nX/cQrQz+zHXsZSkJ+5+g6/d3g7svMEGFbaRwh4Mjr6Kf+Bh9Xa38r3Erh8wCij25VJipXBpVz4++\n8hlcrom1iYrQFeVy8LObL6LDiqLUM4n5x5XiFa5mGx8h/vif6f71GnRrpf+DCttI4Q5l/V3eEZEb\nbqPJHcUT6jZqzDQOe9Ip9WSz0NXMv33pGtLSpr40mAgtH52fyVVzk9jvyabdipr4DpRiq1rD77kJ\no/ME7v9YhbtkozSdTBNSuEOR1nDsNfR/Xo7+4He8a1zGf5m30OyJ46Ang53uGcx0dvCvt6xk5syZ\ndqcVAfLTmy8mxmXw+sA8uqyISe3juJrNY9xBs07E9cLf0PvEp6BVRlmGOyncoebUPnj60/D7z9PR\n2cVvnbezWa/GjZMd7pnsds+gwNnBv1w/hyWLp8/6i+JcmQlRbPjqGkzl5PWBogl3ERzSrhJ5XN/K\na1yJcXI35q8vYWDT96C7xc+JRbAoHYC3TitXrtR79+71+36nLa2hfAvs+i84/hp9Rizb1Cp2mYuw\nlIPTlovt7lnUWwmsiGjgZ3dcTtG8eXanFkGyu6KJL/z3TuJUPx+POEq08kx6X3H6NB9X21msD2Ma\nkeiL7sJ5xTcgcfJrggr/UEoVa619mqtCCredWivg0Ivofb9HNR+j14hjD0t5T19EHxF4tOKAJ4v9\nnixAcVVcHT++96aQWBhBBNfWQ6f4m/8txonFZa4aZjlamcqCQGm6hctVMUv0QRTgnnE5ESu/CPM/\nCRHnX3VeBIbfC7dS6jrgEcABPK61/tn5tpfCPYaBHji5Gyq3o4+/hqovBeCUkUuxsYISz2xM5aRP\nOyk3UzjkyeS0jmSWs50vLo3j1k9fKwNsLmD7q5t44OkdnOhxMsNoY6XrJIlG/5T2mag7uVgdZIk+\nSBKd3rPwmatxFn0cZn8U0uZBGMx7Mx34tXArpRzAMeBjwElgD3Cb1vrQWI+54Au3ZUFXHbRVQuNh\nqC9F15dC/QGUNYCFwSkjhyNqHgetObQTz2krgnorjhNWEjVmEhYG6UY3H8vs4eu3XEdW1iS6hYlp\nx7Q0P3l+J08VN2NikGl0Mc/RTL6jnUhlTn7HWjOTkywyKphtVZBKm/d4rjjIWYEjfyWkz/cW8rS5\nEDmxwUFifP4u3KuAH2itPz749XcAtNY/HesxYV+4tQbTDeaA9+buBU8fuHu8Q877u7y3vnboaYXe\nNu+oxq46dGcddJxAmR+eCfWpGGpJ5wQ5HGEOpWY+rTqBTh1JpxVFq46mW3tX+Y5SHgqdbXxifhI3\nrbuEnJycablIrpia6sY2fv2X3bxR3j28Qnyi6iXd6CZJ9RGrBohVA0QpNxHKJAITh/K9WTRJdzBL\nnSDfaCLLrCWDZhxYw983IxIgMQ8jKR8Vlw6x6RCTBlGJg7cEiIgDVwxExIAzGpwR4IwCR4TMYDiK\niRRuX6Z1zQVOnPH1SeAjkwk2nrJ/vginde5bP8XYT7iRJU0Pb6uGH6eHP1dYKMDAQmmNgcbAwoE5\n/NE5eJuIXh1BM4k06SQadTKndCFVOotqnUGZlUstaeckBXBgkWD0k2H0MCv+NKvmpHH5kiLmzJ5N\nRMTkun+JC8PMjGR+fs/HMU2Tv7x3kNf3VXKkyeRkfyJljN6vX2HhQOPAwkCjFN6PZ7xmvB/Pfb05\n8ZCvGpml6pmpGsg2W8jqbSWz4RApqotkunBN4IzfQuHBgTl8MwZTGZgYaDX0alWDCT9Mqc9IqEe8\nrj78/MyfQI/y2vPdxB7b7Uxg0Xd3TOF4vvHbfNxKqXuBewFmzJgxqX20xsyiu330Lkq+/vLP/aN6\nv7ZGfG5gDZdyY/hJ4xl8ApkYuHEO3/q1iwFc9OOilyh6iKRHR9FFDF3E4MaFQuNQCofDwFDgVOA0\nIN9QFDo7iXYZxEU4SI2PJDMpluykWGZnp5CelnZBzj8t/MPhcHDT5Uu56fKlAPT391NRU0t5fRvV\nTZ3UtZ6ms3eArn4PfW6LAVPjtjQey/vG0sK7XBx4vx5Z8EbqJJUSUinhw26o2gJtWYAmjl7i6CVW\neT9G00+UGiCKfiLwDN7c3pMj5T1B8pZtC8fgK/LD8m2dUbIZ/nj2ydmQsU7uznfSN57JvNHtMaOY\n09dHVNQkBk1NgDSVCCFECJhIU4kvl4v3AHOVUrOUUhHArcBLUwkohBBi8sZtKtFae5RSXwNew9sd\n8Emt9cGAJxNCCDEqn9q4tdavAq8GOIsQQggfSM96IYQIM1K4hRAizEjhFkKIMCOFWwghwowUbiGE\nCDMBmdZVKdUEVPt9x5OTBjTbHWIcoZ4x1PNB6GcM9XwgGf1hKvlmaq3TfdkwIIU7lCil9vo6Gsku\noZ4x1PNB6GcM9XwgGf0hWPmkqUQIIcKMFG4hhAgzF0LhfszuAD4I9Yyhng9CP2Oo5wPJ6A9ByTft\n27iFEGK6uRDOuIUQYlq5oAq3UupvlVJaKTX6EiE2UUr9m1LqiFJqv1LqBaVUkt2ZhiilrlNKHVVK\nlSml/sHuPGdSSuUrpbYqpQ4ppQ4qpb5hd6axKKUcSqkPlFIv251lNEqpJKXUHwefh4cH5+EPGUqp\nbw7+jQ8opf6glArsSgW+ZXpSKdWolDpwxn0pSqk3lFLHBz8mB+LYF0zhVkrlA9cCNXZnGcUbwGKt\n9VK8CzN/x+Y8wPBC0b8BPgEsBG5TSi20N9UIHuBvtdYLgcuAB0Is35m+ARy2O8R5PAJs0lrPB5YR\nQlmVUrnA14GVWuvFeKeXvtXeVAA8BVx31n3/ALyptZ4LvDn4td9dMIUb+CXwd5y7IpPttNava609\ng1/uBPLszHOGS4EyrXWF1noA2ADcaHOmYVrrOq31+4Ofd+EtNrn2pjqXUioP+CTwuN1ZRqOUSgTW\nAk8AaK0HtNbt9qY6hxOIVko5gRjglM150FpvA1rPuvtG4OnBz58GbgrEsS+Iwq2UuhGo1VqX2J3F\nB18G/mp3iEGjLRQdcoURQClVAKwAdtmbZFQP4z1psMbb0CazgCbgfwabcx5XSsXaHWqI1roW+Dne\nd8t1QIfW+nV7U40pU2tdN/h5PZAZiINMm8KtlNo82P519u1G4B+B/xfC+Ya2+S7et//P2Jc0/Cil\n4oDngQe11p125zmTUupTQKPWutjuLOfhBC4CHtVarwC6CdBb/MkYbCe+Ee8/mBwgVil1h72pxqe9\nXfYC8g7fb6u8201rfc1o9yulluD9g5co77LNecD7SqlLtdb1ducbopS6G/gUcLUOnT6atUD+GV/n\nDd4XMpRSLrxF+xmt9Z/szjOKNcANSqnrgSggQSn1O611KBWek8BJrfXQu5U/EkKFG7gGqNRaNwEo\npf4ErAZ+Z2uq0TUopbK11nVKqWygMRAHmTZn3GPRWpdqrTO01gVa6wK8T9KLglm0x6OUug7vW+kb\ntNY9duc5Q0gvFK28/4mfAA5rrX9hd57RaK2/o7XOG3zu3QpsCbGizeBr4YRSqmjwrquBQzZGOlsN\ncJlSKmbwb341IXTx9CwvAXcNfn4X8GIgDjJtzrjD3H8AkcAbg+8Kdmqtv2pvpLBYKHoNcCdQqpTa\nN3jfPw6ukSom5v8Czwz+g64AvmRznmFa611KqT8C7+NtSvyAEBhBqZT6A7AOSFNKnQS+D/wMeE4p\ndQ/eGVI/H5Bjh867ciGEEL6Y9k0lQggx3UjhFkKIMCOFWwghwowUbiGECDNSuIUQIsxI4RZCiDAj\nhVsIIcKMFG4hhAgz/z+i8FHoqCRn6AAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f3858207fd0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "def func(x):\n",
    "    return np.exp(-x ** 2) + 3 * np.exp(-(x - 3) ** 2)\n",
    "x = np.linspace(-5, 10, 100)\n",
    "#提议分布\n",
    "rv = Gaussian(mu=np.array([2.]), var=np.array([2.]))\n",
    "plt.plot(x, func(x), label=r\"$\\tilde{p}(z)$\")\n",
    "plt.plot(x, 15 * rv.pdf(x), label=r\"$kq(z)$\")\n",
    "plt.fill_between(x, func(x), 15 * rv.pdf(x), color=\"gray\")\n",
    "plt.legend(fontsize=15)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAD8CAYAAACMwORRAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xd8m+W99/HPT5Yt79iOZzyzByHTLSEUyigcdqDQEkbK\namkKdDysFnoKpX066GlPD7Q80DBaWlKgrMKBMAKEGRISJziD7OEVr3jb8pSu5w9Jju3YsWJLliz/\n3q+XXrGkW/f988hXl677uq9LjDEopZQKLZZAF6CUUsr3NNyVUioEabgrpVQI0nBXSqkQpOGulFIh\nSMNdKaVCkIa7UkqFIA13pZQKQRruSikVgqyBOnBycrLJy8sL1OGVUmpUKigoOGyMSRlsu4CFe15e\nHhs3bgzU4ZVSalQSkSJvttNuGaWUCkEa7kopFYI03JVSKgRpuCulVAgaNNxFJFJEPhORQhHZLiL3\n97ONiMhDIrJXRLaIyAL/lKuUUsob3oyWaQfONMY0i0g48LGIvGGMWddjm/OAqe7bScAj7n+VUkoF\nwKDhblxLNTW774a7b32Xb1oC/N297ToRSRCRDGNMuU+rVUr1q7GxkaqqKjo7OwNdihqG8PBwUlNT\niY+PH/a+vBrnLiJhQAEwBXjYGLO+zyaZQEmP+6XuxzTcVcgqb2jlk701fH1+JhaLBKyOxsZGKisr\nyczMJCoqCpHA1aKGzhhDa2srZWVlAMMOeK/C3RjjAOaJSALwsojMNsZsO96DichNwE0AOTk5x/ty\npYJGYUk93/77Rqqb2omOCOP8EzMCVktVVRWZmZlER0cHrAY1fCJCdHQ0mZmZHDp0aNjhflyjZYwx\n9cAa4Nw+T5UB2T3uZ7kf6/v6FcaYfGNMfkrKoFfPKhWUVm0t55t/+RSb1ULe+Gj+uHo3DmfgFprv\n7OwkKioqYMdXvhUVFeWT7jVvRsukuFvsiEgUcDaws89mrwLfco+aWQQ0aH+7CkVvba/g5pWbmJ05\njn/fcgp3/Md09lQ189qWQwGtS7tiQoevfpfetNwzgDUisgXYAKw2xrwmIstFZLl7m1XAfmAv8Bhw\ns0+qUyrIvFhQSsa4SFZ++ySSY22cPzuDGelxPPjOHroczkCXp1Q3b0bLbAHm9/P4oz2+NsAtvi1N\nqeDS0eVk7b4aLpo7gcjwMAAsFuFHX5vG8qcLeLXwEF9fkBXgKpVy0StUlfLSpuI6mtu7OH167/NF\n/3FCGidMiOfBd/fQqa13FSQ03JXy0vu7qrFahMWTx/d6XES4+fQpFNXY2XiwLkDVhY5rr72WxMRE\nHnzwQZ555hmuu+66Xs8bY5g3bx5PPfWU1/u89dZbufHGG31caXDTcFfKSx/srmZhbiJxkeFHPecJ\n/E3FGu7DUVdXx1tvvcXjjz/O008/zS233MINN9zQa5t//etf1NbWctVVV3m93zvuuIOVK1eyd+9e\nX5cctDTclfJCZWMbO8ob+er0/ofwJsZEMCk5hs0a7sOSmJhIRUUFl112GRs2bKC2tpbTTjut1zYP\nPfQQy5YtIzz86DfZgeTl5fGVr3yFRx55xNclBy0Nd6W88MHuagBOn5Y64DbzcxLZXFyPa3yBGoqk\npCR+9rOfcfXVVzN+/HiSkpL47ne/S3t7OwB79+5l7dq1XH755b1e99xzzyEi/d4aGxsBuOyyy1i5\nciVO59g4L6LhrpQXPthdTUqcjZkZcQNusyA3gZqWDopr7SNYWegoLi6mrq6OP/7xjyQlJfH888/z\nox/9iBUrVvDwww8D8O677xITE8PcuXN7vfa0007j008/7b69/vrrJCUlcf7553df6bl48WIqKyvZ\nunXriH9vgRCwNVSVGi26HE4+3nOYs2elHfMCkwU5iYCr3z13fMxIlRcyCgsLAdfJz9/+9rcAnHnm\nmXz44Ye8+eab3HbbbRQUFDBz5kwslt7t0oyMDDIyXFNANDY2ctZZZzF79myef/757m1OOOEEwsLC\n+Oyzz456cwhFGu5KDaKwtIGG1k6+Ou3YU2ZMS4sjJiKMzcX1XDo/sOPd7//f7XxxqDEgx541IZ77\nLjrhuF9XWFhIdHQ0d955Z6/Hp02bxmeffQZARUUFycnJA+6jpaWF888/H4vFwmuvvdZrvh2r1UpC\nQgIVFRXHXdtopOGu1CA+2lONReDUqQOHCkCYRZibnaAjZoaosLCQ0047jfHjew81raioYMKECQC0\ntbUNOEFaW1sbF198MS0tLbz33nvExR3dhWaz2Whra/N98UFIw12pQWw/1MjE5BgSoiMG3XZBTiKP\nfLAPe0cX0RGB++81lJZzoBUWFnLqqaf2eqypqYn33nuP++93LQCXlJTUb8u7s7OTyy67jPLycj74\n4AMSExP7PUZ9fT1JSUm+Lz4I6QlVpQaxu7KJ6ekDn0jtaUFuAg6nYWtpg5+rCi0tLS3s27ePffv2\n9Xr8oYceAlwXNgFMnz6dAwcO9NrG4XBw5ZVXsmvXLt555x0GmnG2uroau93OtGnT/PAdBB8Nd6WO\nwd7RRXGtnelp3s2tPS/bc1K13p9lhZytW7fidDrZt28fv/jFL1izZg133303999/PytWrCAhIQGA\nU045heLiYqqrq7tf+73vfY9Vq1bx61//muLiYtatW8e6deuOGhWzceNGRITFixeP6PcWKNoto9Qx\n7KlsxhiYnh7r1fZJMRFMTI7RfvfjVFhYSHx8PK+//jrXX389v/rVr5g+fTrPP/88S5Ys6d7u9NNP\nJykpiTfffJNly5ZhjOHZZ5+ltbWVK664otc+L7roIl599dXu+2+++SZf/epXj+rTD1Ua7kodw67K\nJgCmp3u/Ks78nAQ+3H0YY4zOs+6lwsJCTjzxRObMmUNBQcGA20VERHDNNdfw7LPPsmzZsl4XKR2L\nw+HgxRdf7B5iORZot4xSx7Crogmb1UJOkvdL2C3ISeRwczulda1+rCy0FBYWej32/M4772TNmjXs\n3r3b6/0///zzREVFsXTp0qGWOOpouCt1DLsrm5iaFkvYcSyAPSdrHADbD+lJVW8YY9i6dStz5szx\navusrCyefPJJysu9X+zNGMMTTzyB1Tp2OivGzneq1BDsrGjitKnHt97v5BRX//y+6hZ/lBRyvO1a\n6el4W+BXXnnlcW0fCrTlrtQAals6qG5q9/pkqkeMzUp6fCT7qpv9VJlSg9NwV2oAu4dwMtVjUkoM\n+7XlrgJIw12pAeyqcId7mncXMPU0OSWW/dXNOv2vChgNd6UGsKuyiXFR4aTF2477tZNSYmhs6+Jw\nc4cfKlNqcBruSg1gV0UT09PihjRW3XNSdb/2u6sA0XBXqh/GGHZXNDHtOE+mekxKcc3nriNmVKAM\nGu4iki0ia0TkCxHZLiI/7Geb00WkQUQ+d9/u9U+5So2M8oY2mtq7hnQyFWDCuCgiwy3aclcB4804\n9y7gdmPMJhGJAwpEZLUx5os+231kjLnQ9yUqNfKGczIVwGIRJibH6nBIFTCDttyNMeXGmE3ur5uA\nHUCmvwtTKpC655QZYrgDTE6JYf9h7ZZRgXFcfe4ikgfMB9b38/RiEdkiIm+IyOhbKUCpHnZXNpEW\nb2NcdPiQ9zEpJZaSWjvtXQ4fVqaUd7wOdxGJBV4EfmSM6Xut8CYgxxgzB/gT8O8B9nGTiGwUkY09\n52NWKtgU1djJG+Yi15NTYnAa177U6Pf+++8jImzbti3QpXjFq3AXkXBcwb7SGPNS3+eNMY3GmGb3\n16uAcBE5asFJY8wKY0y+MSZ/oNVSlAoGRTV2csd7PxNkf7rnmKnSfnc18rwZLSPAE8AOY8x/D7BN\nuns7ROTL7v3W+LJQpUZKS3sXh5vbyR1my31isuv12u+uAsGblvspwDLgzB5DHc8XkeUisty9zeXA\nNhEpBB4Clhq97lqNUsW1rm6U4bbcY2xWMsZFasvdS9u3b+fcc88lKSmJmJgYZs6cycMPPwzA66+/\nztlnn01qairx8fEsWrSIt99+u9frf/7zn5OcnMz69evJz88nKiqKr3zlKxw4cICqqiouueQSYmNj\nmTlzJu+9916v1+bl5XHHHXfwy1/+kvT0dGJjY7n66qtpaDj2tM1Op5Pf/va3TJkyBZvNxrRp03jq\nqad6bfPxxx9z6qmnEh8fT3x8PPPmzeP555/3wU/s2AYdCmmM+Rg45iV6xpg/A3/2VVFKBZKnjzw3\naXgtd3BdzLRPW+5eueiii5g5cyZPP/00NpuNXbt2dU8FfODAAc4//3xuv/12wsLCeOONNzjvvPP4\n8MMPOeWUU7r3Ybfbuemmm7jrrruIiYnhBz/4AcuWLcNms3Heeedx880387vf/Y5vfOMblJSUEB19\n5A38mWeeYcqUKTz22GOUl5dz11138e1vf/uYQfz973+fp556invvvZcFCxawevVqbrjhBsaPH8+F\nF15IY2MjF154IUuWLOHee+/tnru+vn4E1tg1xgTktnDhQqNUMHr0/b0m98evmXp7x7D39bN/bzWz\n733TOJ1OH1TWvy+++MKn+3M4nKaqsc2vNfdVXV1tALNly5ZBt3U4HKazs9Occ8455vrrr+9+/L77\n7jOAef/997sfe/jhhw1g7r///u7Htm/fbgCzatWq7sdyc3NNYmKiaWpq6n7s6aefNiLS/fNds2aN\nAczWrVuNMcbs2bPHiIj529/+1qu+ZcuWmfz8fGOMMRs2bDCAaWxsPJ4fxzF/p8BG40XG6vQDSvVR\nVGsnITqccVFDHwbpMSk5hqb2Lqqb231Qmf85nYYrH1vHyb95l6Ur1uF0jkzvalJSEtnZ2Sxfvpzn\nnnuOqqqqXs+XlpZy7bXXkpmZidVqJTw8nLfffvuopfYiIiI49dRTu+9PmTIFgDPPPPOox8rKynq9\n9uyzzyY29sh0E5deeinGGDZs2NBvze+++y4Wi4VLL72Urq6u7ttZZ53F559/jsPhYPLkycTGxnLV\nVVfxyiuvjEyL3U3DXak+imvswz6Z6jE51TNiZnR0zdS0dFBQVEeX01BQVEdNy8jMammxWHj77bdJ\nT0/nhhtuID09nVNPPZXNmzfjdDq5+OKLWbt2Lb/4xS9Ys2YNGzZs4LzzzqOtra3XfuLi4rBYjsRa\nREQEAAkJCUc91ve1qampve5HR0cTGxs74HJ+hw8fxuFwMG7cOMLDw7tv1113HV1dXZSXl5OYmMjq\n1avp7Ozkm9/8JikpKVxwwQXs379/6D8sL+kye0r1UVTbwvzsRJ/syzNWvqimhZMnj/fJPv0pOTaC\nhbmJFBTVsTA3keTYiBE79owZM3jxxRfp7Ozko48+4sc//jEXXHAB77//Pps3b+aNN97g3HPP7d6+\ntdW3C5D3/bRgt9tpbm4mIyOj3+2TkpKwWq188sknvd5QPDxvFosWLeLNN9+ktbWVd955h9tuu42r\nrrqKdevW+bT+vjTcleqho8tJWV0rl8zzzQwbGeMisVqEkrrRcSGTiPDMdxZR09JBcmzEkKY7Hq7w\n8HDOPPPM7hD0tJxttiPz6hcVFfHJJ594vai2N1avXk1zc3N318zLL7+MiJCfn9/v9meeeSYOh4OG\nhgbOPvvsQfcfFRXFRRddxLZt2/jNb37js7oHouGuVA9l9a04DeQkDW8YpIc1zMKEhCiKa33byvQn\ni0VIiTv+BUqGY8uWLdxxxx1cccUVTJo0ibq6Oh544AHmzp3LokWLyMrK4vbbb+eXv/wlTU1N3Hff\nfWRm+naKq6ioKC644ALuvPNOysvLufPOO7n00kuZNWtWv9tPnz6d5cuXs3TpUu666y7y8/Npa2tj\n+/bt7N69m8cff5zXX3+dJ598kksuuYScnBzKysr4y1/+0uscgL9ouCvVQ1GNq288L9k3fe4A2UlR\nlNSOjpZ7oKSnp5OWlsavfvUrDh06REJCAmeccQYPPPAANpuNl156iVtuuYXLL7+crKwsfvrTn/L+\n++/7dCqApUuXEhcXx4033khzczMXX3wxjzzyyDFf8/DDDzNt2jQee+wx7r33XuLj45k1axY33ngj\n4Dp5KyLcc889VFVVkZKSwoUXXsivf/1rn9U9EDEButYoPz/fbNy4MSDHVmogf//0IPe+sp3P7jmL\n1PhIn+zz7pe2sPqLSjb+5+Af3Ydix44dzJw50y/7Hivy8vK4/PLL+f3vfx/oUoBj/05FpMAY039f\nUQ86WkapHopq7ESFh/m0WyIrMZrDzR20tHf5bJ9KDUbDXakeimrs5CRF+/REYra7/760bvT0u6vR\nT/vcleqhqKale8IvX/GcnC2utTM9feiLfyj/OXjwYKBL8DltuSvl5nQaimuHP9VvX9mJUQB6UlWN\nKA13pdyqmtpp73KS46OrUz2SYiKIiQjz61j3QA2MUL7nq9+lhrtSbp5hkLk+GuPuISJkJ0X7reUe\nHh7u86s1VeC0trYSHj78eY003JVyK3KH73CX1+tPVmI0JX66kCk1NZWysjLsdru24EcxYwx2u52y\nsrKj5rkZCj2hqpRbUU0LVoswIcE349t7ykmKZu2+wxhjfH5Jf3x8PACHDh2is7PTp/tWIys8PJy0\ntLTu3+lwaLgr5VZUYyczMQprmO8/0GYnRWHvcLjnbPH9pf2eVX6U8tBuGaXcSutafTanTF/Zia79\n6ogZNVI03JVyK62zk+UetuhrOeOPjHVXaiRouCsFtHY4ONzcQVaif1runjcNvUpVjRQNd6WAsnpX\ni9pfLffoCCvJsRHaLaNGjIa7UtA9TNFfLXdwzTGj3TJqpGi4K4Wrvx2OTBXgD9mJ0aNmRSY1+g0a\n7iKSLSJrROQLEdkuIj/sZxsRkYdEZK+IbBGRBf4pVyn/KK1rJcJq8cswRY+cpGgO1bfR5XD67RhK\neXjTcu8CbjfGzAIWAbeISN91p84DprpvNwHHXr5EqSBTWtdKVmIUFov/1gzNTorC4TSUN7T57RhK\neQwa7saYcmPMJvfXTcAOoO/ihUuAvxuXdUCCiPS/ZLhSQaikzu7X/nbQse5qZB1Xn7uI5AHzgfV9\nnsoESnrcL+XoNwClgpan5e5P2Uk61l2NHK/DXURigReBHxljGodyMBG5SUQ2isjG6urqoexCKZ9r\nae+itqXD7+GeMS6SMItQVq9j3ZX/eRXuIhKOK9hXGmNe6meTMiC7x/0s92O9GGNWGGPyjTH5KSkp\nQ6lXKZ/zhK2/u2WsYRbS4yP1QiY1IrwZLSPAE8AOY8x/D7DZq8C33KNmFgENxphyH9aplN94+sD9\nOQzSIysxqnvYpVL+5M2skKcAy4CtIvK5+7F7gBwAY8yjwCrgfGAvYAeu932pSvmHpyXt75a75xhr\n9x32+3GUGjTcjTEfA8ccH2ZcKwTc4quilBpJpXV2bFYLybERfj9WZmIUFY1tdHQ5ibDqNYTKf/Sv\nS415JbWukTK+XkSjP1mJURgD5Q3a7678S8NdjXml9fbuYYr+5hmRU6YnVZWfabirMW8kxrh7eC5k\n0hEzyt803NWY1tTWSb29c0ROpgKkj4vEIuiIGeV3Gu5qTDsyUmZkWu7hYRYyxkVpy135nYa7GtM8\nIZs9Qi13gMwEDXflfxruakzzdI+MVMvdcyydgkD5m4a7GtNK61qJCg8jKcb/Y9w9shKjKG9opVPn\ndVd+pOGuxrSSWvuIjXH3yEqMxmmgQud1V36k4a7GtJIRHAbp4TmeLrmn/EnDXY1ppXUjdwGTR5aO\ndVcjQMNdjVkN9k6a2rpGdKQMuMa6i2i4K//ScFdjlqdbJDtpZLtlIqyued11CgLlTxruaszyzOM+\nUlen9qTzuit/03BXY9aRlnsgwj1au2WUX2m4qzGrpLaV+Egr46LCR/zYWe553bt0rLvyEw13NWaV\n1tkD0iUDrikIHE5DuY51V36i4a7GrJK61hE/merheVPRaQiUv2i4qzHJGOMa4x6glrvnQibtd1f+\nouGuxqTq5nbaOp0BOZkKkJHgGuvuGbGjlK9puKsxqaTWPdVvgLplbNYw0uMjdQoC5Tca7mpM8owx\nD1S3jOfYpbXaLaP8Q8NdjUlHVmAKXLhnJUVpy135zaDhLiJPikiViGwb4PnTRaRBRD533+71fZlK\n+VZJrZ3kWBtREWEBqyEnKZqKxjbauxwBq0GFLm9a7n8Dzh1km4+MMfPct18Mvyyl/Kukzj7iU/32\nlZ0YjTHoHDPKLwYNd2PMh0DtCNSi1IgpqW0N2EgZD8/xSzTclR/4qs99sYhsEZE3ROQEH+1TKb9w\nOA2H6lvJDnTL3T1Sp1iHQyo/sPpgH5uAHGNMs4icD/wbmNrfhiJyE3ATQE5Ojg8OrdTxK29opctp\nAt5yT4uLJCLMQqmGu/KDYbfcjTGNxphm99ergHARSR5g2xXGmHxjTH5KSspwD63UkHhGygRyGCSA\nxSJkJeqIGeUfww53EUkX9+rCIvJl9z5rhrtfpfzFc1VooC5g6ikrKVq7ZZRfDNotIyLPAKcDySJS\nCtwHhAMYYx4FLge+JyJdQCuw1Bhj/FaxUsNUUteKRSBjXODDPScpisKS+kCXoULQoOFujLlykOf/\nDPzZZxUp5WeltXbS4yOJsAb+Gr7sxGgaWjtpbOskPnLk55VXoSvwf91KjbCSOjtZAT6Z6tE9HFK7\nZpSPabirMae41k5OsIR7ooa78g8NdzWm2Du6qGxsZ2JyTKBLAeh+kynRCcSUj2m4qzGlqMbVQs4d\nHxwt93HR4cRFWnU4pPI5DXc1phw83AJA3vjgaLmDq2tGh0MqX9NwV2PKQXfLPS9IumXANd5e+9yV\nr2m4qzHl4OEWkmNtxNp8MfOGb+QkRVNa14peHqJ8ScNdjSkHa1rIC5L+do/spGjau5xUN7UHuhQV\nQjTc1ZhysKYlqLpk4MhwSO13V76k4a7GDM8wyGBsuQM6Ykb5lIa7GjOKgvBkKtC9IlRxjY51V76j\n4a7GjKKa4BsGCRAZHkZ6fCRFtS2BLkWFEA13NWYcOBxcFzD1lJcc3T0GXylf0HBXY4ZrGGQEcUE4\n++LE5FgOaLgrH9JwV2PGwZoWcoOsS8ZjYnI0dfZO6u0dgS5FhQgNdzVmuMa4B2u4xwJo6135jIa7\nGhOCdRikh2eWSg135Ssa7mpMCNZhkB45SdFYBD2pqnxGw12NCcE6DNIjwmohKzGa/Rruykc03NWY\n0D0MMjk4u2XA1TVzsEbDXfmGhrsaE4pqWhgfExHUi1BPTI7hQHWLzg6pfELDXY0JBw4H34RhfU1M\njqGlw6GzQyqf0HBXY4JrjHvwdsnAkZO9OmJG+cKg4S4iT4pIlYhsG+B5EZGHRGSviGwRkQW+L1Op\noWuwd1LZ2M7U1LhAl3JMkzTclQ9503L/G3DuMZ4/D5jqvt0EPDL8spTynd1VTQBMT48NcCXHNiEh\niogwCwf0pKrygUHD3RjzIVB7jE2WAH83LuuABBHJ8FWBKvgZY9hSWk9xTXDOR76rwhPu8QGu5NjC\nLELO+GgOVGu4q+HzxUKSmUBJj/ul7sfK+24oIjfhat2Tk5Pjg0OrQDLG8N7OKh5es5dNxfWIwLkn\npPPdr05mXnZCoMvrtruyiViblQnjIgNdyqAmJsdot4zyiRFdJdgYswJYAZCfn6/jvUax1g4H1zyx\nnoKiOjITovj5RbOoamrnH+uKeGNbBUu/lM1vvn4iIhLoUtlV0cS0tNigqGUwE5Nj+GB3NQ6nIcwS\n/PWq4OWLcC8Dsnvcz3I/pkLY/f+7nU3Fdfzq0tl8Mz+b8DBXD9/NZ0zhD2/v4q+fHGROVgJXnRTY\nT2jGGHZVNnHe7PSA1uGtickxdHQ5OVTf2r38nlJD4YuhkK8C33KPmlkENBhjjuqSUaHjlc/LeHZD\nCTefPpmrT8rtDnaAWJuVn10wi9OmpfDz/93OtrKGAFYK1U3t1Ns7mZ4W3CNlPDwTiOmVqmq4vBkK\n+QzwKTBdREpF5EYRWS4iy92brAL2A3uBx4Cb/VatCrgDh1u456Wt5Ocm8n++Nq3fbSwW4Y/fnEtS\ndAS3/nMTjW2dI1zlEbsqXSdTp6WPrnDXfnc1XIN2yxhjrhzkeQPc4rOKVNByOA3ff2YT4VYLD105\nH2vYwG2D8bE2/nTVfJauWMd9r2znj1fMG8FKj+geKTNKWu6pcTaiI8LYryNm1DDpFarKa69tOcS2\nskZ+sWQ2ExKiBt3+S3lJfPe0Sby8uaw7ZEfa7somkmMjGB9rC8jxj5eIMDkllr1VzYEuRY1yGu7K\nK06n4eE1e5mWFsuFJ3p/GcNNp00i1mbloXf3+LG6ge2qbGb6KOmS8ZiRHsfOisZAl6FGOQ135ZW3\ntlewu7KZW86YguU4huglREdw/Sl5vL61fMQDy+k07KlsYtoo6ZLxmJ4ex+HmDg436wRiaug03NWg\njDH86b29TEqO4cI5E4779Td+ZSJxAWi9l9a1Yu9wjJr+do+ZGa4raQPVlaVCg4a7GtS7O6r4oryR\nm8+YMqQLazyt91VbK0a09T7aRsp4eLqRdpRr14waOg13dUyuVvsespOiWDLv+FvtHjcEoPW+2xPu\no6zlnhxrIznWpi13NSwa7uqYCorqKCxtYPlXJ/e6WOl4JURHcM3Juby5rYKy+lYfVjiwXRVNZCVG\nEWsb0Vk2fMJ1UlXDXQ2dhrs6pn+uLybOZuXS+ZnD3tc1i3IBWLmuaNj78sbuyqZR19/uMSM9jt2V\nTTicOgWTGhoNdzWgensHr20t55L5mURHDL/1m5kQxVkz03huQwntXQ4fVDiwtk4H+6qbR11/u8f0\n9Djau5w6DYEaMg13NaCXNpXR0eXkyi/7bvKvb52cS01LB6u2+nf6oa1lDXQ6DAtyEv16HH/RETNq\nuDTcVb+MMTzzWTFzsxOYNcF3i1ycMjmZSckx/ONT/3bNbDxYB8CCnOCZV/54TEmNxSKwU0fMqCHS\ncFf9KiiqY09VM1f7sNUOrknFrlmUy6bier/OGFlQVMuk5JhRM+1AX5HhYUxMjtGTqmrINNxVvzwn\nUi+c6/sVEy9bmEVUeJjfWu/GGAqK6liYOzq7ZDxmpMdruKsh03BXR2mwd/K6D0+k9jUuKpxL5k/g\nlcIyGlpYH5exAAAYi0lEQVR9Px3wvuoW6uyd5OeN9nCPo7jWTkt7V6BLUaOQhrs6yquFZbR3Obni\nS9mDbzxEV5+US1unk5c3lfp83wVFrvXcF+Ym+XzfI8lzparnSluljoeGuzrKC5vKmJEex+zMcX47\nxuzMcczNGsfK9cW4lgTwnY0H60iMDmdySoxP9zvSZqTriBk1dBruqpe9VU0UltRz+cIsvx/r6pNy\n2VPVzMaiOp/u19PfPhoWxD6WrMQoYiLCdMSMGhINd9XLCwVlhFmEJfOGf0XqYC6cm0GczerTK1Zr\nmtvZf7hl1HfJgGtk0YyMeLYGeB1aNTppuKtuDqfh5c2lnDE9hZQ4/w8hjI6w8vUFmazaVkFtS4dP\n9lng/hQw2k+meizMTWRbWaPfr+hVoUfDXXX7eO9hKhvbuWyB/7tkPK46KZeOLicvFvjmxGpBUR0R\nYRZO9OP5gpG0ICeRDofTr9cEqNCk4a66vVBQSkJ0OGfOTB2xY05PjyM/N5F/flaM0weTZG0sqmN2\nZjyR4WE+qC7wPGP1PVfcKuUtDXcFQGNbJ29vr+DiuROwWUc2GJednMuBwy2s2VU1rP20dTrYWtpA\nft7o72/3SImzkTs+uru7SSlvabgrAP638BDtXc4R7ZLxOP/EDDLGRfL4RweGtZ/3dlbR4XBy2tQU\nH1UWHBbmJlJQVOfzIaMqtHkV7iJyrojsEpG9IvKTfp4/XUQaRORz9+1e35eq/OlfG0qYkR7HnKyR\n76sOD7Nw/Sl5fLq/Zlh9y698XkZKnI2TJ4/3YXWBtzA3kZqWDopq7IEuRY0ig4a7iIQBDwPnAbOA\nK0VkVj+bfmSMmee+/cLHdSo/2lnRSGFpA9/Izw7Y2PClX84h1mbl8Y/2D+n1Da2drNlZzUVzJgxp\nnddglu8e1unr6wFUaPOm5f5lYK8xZr8xpgN4Flji37LUSHp+YynhYeKT1ZaGKj4ynCu+lM1rW8op\nbzj+Zfje2lZBh8PJxcNY5zVYTU2NJS7Sqv3u6rh4E+6ZQEmP+6Xux/paLCJbROQNETmhvx2JyE0i\nslFENlZXVw+hXOVrHV1OXt5cxtmz0kiKiQhoLdefkocB/vbJweN+7auFh8gdH83cAHQr+ZvFIizI\nSeyeM0cpb/jqhOomIMcYMwf4E/Dv/jYyxqwwxuQbY/JTUkLrpNdo9e6OSmpbOvhGvv8mCfNWVmI0\n581O55/ri6k7jouaqhrbWLvvMEvmThj1Uw4MJD83kd2VzX6ZRVOFJm/CvQzo+T8/y/1YN2NMozGm\n2f31KiBcRJJ9VqXym+c2lpAeHxk0I0y+f+ZU7J0O/rB6l9eveW1LOU5DSHbJeHjGu28q1q4Z5R1v\nwn0DMFVEJopIBLAUeLXnBiKSLu4mk4h82b3fGl8Xq3yrvKGVD3dXc/nCrKA5CTk9PY5vnZzLyvXF\nXo+ceaXwELMy4pmSOjoXw/bG3OwEwizCJu13V14aNNyNMV3ArcBbwA7gX8aY7SKyXESWuze7HNgm\nIoXAQ8BSo4Nyg96zn5VggG/kj/zY9mP50demMT4mgp+9sm3Qq1YLimopLKlnSQi32gFibFZmZsSx\n/oD2uyvveNXnboxZZYyZZoyZbIz5lfuxR40xj7q//rMx5gRjzFxjzCJjzFp/Fq2Gr6PLycr1xZwx\nPZXc8cE17/m4qHB+ct5MNhfX88IxFvNoae/itn8VkpkQxZUn+Xat12B02tQUCorqqLf7ZpI1Fdr0\nCtUxatXWcg43t3Pt4rxAl9Kvr8/PZEFOAg+8sZPiAS7e+dWqHRTX2vnDN+cSHxk+whWOvHNOSMfh\nNLy3c3jTNKixQcN9jPrr2oNMSonh1CnBed7bYhF+/fUT6XIaLv1/nxx1InHNzir+ub6Y75w6iUWT\nQuuK1IHMyRxHWryNt7dXBroUNQr4fvVjFfQ2F9dRWFLP/RefgCVITqT2Z0Z6PC/dvJjr/7qBK1es\n4/9eMpu4yHC+KG/kn+uLmJEex+3nTDv+HR/a7N12E+YHZn8DsFiEs2el8WJBGW2djpCZ+VL5h7bc\nx6Cn1h4k1mblshFYSm+4JqfE8vLNizlhQjx3vrCF5U8X8Of39pAca+PBpfNHfAbLQDtnVjqtnQ4+\n3nM40KWoIKct9zGmqqmN17eWc/VJucTaRsevf3ysjX9+ZxHv76oiLT6SGenxREWMQKh72yIfQYsm\njSfOZuXtLyr42qy0QJejgtjo+N+tfOavnxyky2n41sm5gS7luESGh3Hu7IxAlxFwEVYLZ8xI5Z0d\nVTicJmiuT1DBR7tlxpDDze08tfYgF82ZwKSU2ECXo4bonBPSqG3p0InE1DFpuI8hf/lgH22dDn74\ntamBLkUNw1enpRARZuHt7RWBLkUFMQ33MaKqsY2/f1rEJfMzmayt9lEtLjKcxVPGs2prOV0OZ6DL\nUUFKw32M+H/v76PLafjhWdpqDwVXn5TLoYY23timrXfVPw33MaC8oZV/flbM5Quygm6qATU0Z81I\nZVJyDCs+3K9rq6p+6WiZMeCBN3ZijOHWM6cEupSxzYcXO1kswo2nTuSnL29j/YHaMXOVrvKettxD\n3Hs7K/n354e49YypZCdFB7oc5UOXLcgiKSZiyOvOqtCm4R7Cmto6ueelbUxPi+N7p08OdDnKxyLD\nw1i2KJd3dlSxt6o50OWoIKPhHsJ++8ZOqpraeODyOURY9VcdipadnIvNauGJj7X1rnrTPvcQtXbf\nYVauL+Y7p05kXnbCyB34eC7ZH+ZEWgqSY218Iz+LZz8r4eqTcpmdGXoLhKuh0XAPQQcOt3DLyk1M\nSo7htrOn+2anQTjPinK545zpvL29ktv+9Tmv3voVnS1SAdotE3ION7dz3V8/wyLCX6//0shMsKUC\nKiE6gt9dPofdlc384W3vFxZXoU1b7iGktcPBt5/aSGVjG898Z9HYG9M+hj9dnD49lWsW5fD4xwc4\na2aaDo1U2nIPFTXuFnthaT0PLp3P/JzEQJekRtg9588kNyma2577nKKalkCXowJMW+4hYGtpA9/9\nx0ZqWjr4nyvm8R8npAe6JO+M0ApGwcLphIa2ThKiwhE/zNQbHWHlz1ct4Jon1nPZI2v52/Vf1hOs\nY5i23EexToeTp9Ye5LJH1yIivPi9xSyZlxnoskYNpxPq7J0MdPV+z+cH23aw13c5DHe9WMh1T67n\n7pe34vTTfF+zM8fxwvKTsVnDuOIvn+qKTWOYVy13ETkXeBAIAx43xvy2z/Pifv58wA5cZ4zZ5ONa\nlZvTafjfLYf479W7Kaqxc+rUZP7ninmMj7Ud/86G0E/dswVqTP9f+6Jl2n0cY5Ah7nCg1rLTCff8\neys7DjUwNS2O31x6Is0dju7tup8vb2RKSizWMGFnRRMzM+L5v0tm09DaiQgkRkd0f9/xkVbq7a7H\nx0WF89N/b2NHeSMz0uPochh2VTYB8EVZAw1tnSRGhx9Vb5fDsLe8kWlpsVgsQ2t7TUmN48XvLea6\nv37Gt55cz+ULs/jR16YxISFqSPtTo9Og4S4iYcDDwNlAKbBBRF41xnzRY7PzgKnu20nAI+5/lY8Y\nY9h+qJFVW8t5fWs5RTV2ZqTH8eR1+ZwxPXXI4TeYvuHYN/TCLLCrookZ6fEgsLO8kalpcTxw2RzC\nLNL9+vhIK/WtnQiuQOxbbs/tGtu6iI+0dodjZI6DZ76zaMDFvJ1OQ01LB8mG/gPcHbB3/ccMkmJc\nx65r7eCLsgacwM6KJq56Yj0dXc7u8C6ps7OjvBGH80goA+w41MCPX9zS/dgJGfGAYUd5E5HhYdg7\nHQBMT4tlb1UzDuP6mTh7tPgjI8KIjzz6v16Xw3D1E+tZ31ZOXKSVzf95NtYhXnyWPi6Sfy0/mQff\n2cM/Pi3i358f4uqTclgyL5M5meOCemF05RvetNy/DOw1xuwHEJFngSVAz3BfAvzduKanWyciCSKS\nYYwp93nFIc7pNBxuaaeioY1D9W3srGhkW1kjW8vqqWxsJ8wiLJ48ntvPmc6FJ2YM6z+p02moa+k4\nZuB6wnFmRjy/vuREGto6+w297eWNCGBwheVdL2zhgcvm8J+vuALaZrVg73AF3wkTxvGbS0/E0zB1\nOuHul7eys6KRSGsYbZ1dTE2LY09VMw6noa24gIb9YQO2dH/80hb2VjV31wiuNwpjTHet2w81csNT\nG5iSGstvLz2R3721i549I63u2naUN/Kj5zZTXGvHZg2j1eno3iZMhCmpMb2+7y/KG/HktifYAfZU\nNTPN/T1MT41hf42dtk7XEds6Ha6fuwjG/WoBSmtbaWnvAqCprYu91c3MyIj34jfZv/jIcH524Sxu\n+MpE/mf1bv7+aRF//eQgybE2vjothVkT4pmSGsuk5BhS421jbrHxUOdNuGcCJT3ul3J0q7y/bTIB\nn4f75uI6/v5pka93OyQ9p1o1vR4/ct8Y4+qzNQanMTic4HA66XIa2rucrlung+b2Lpraumhq6+zV\nyhOBSckxnDxpPCdPHs85s9JJjIkYdu1Op2Hpik+xFxUARwcu0CvIdxxqoNbejkWEKSmxvQLOY1Jy\nNPsO2wHYVdlEcW1L9+s9wQ6ulmx9awdJ7u+jzt7B9kMNALR0uMJtT2UTU9LiukM7IeroYHc64Scv\nbWFXhauWHeWN1LV28F9v7WJHeSN5yTFMT41hV6WrBe1wGnZVNHHHi1s4MMBcLMYYDta4voe2HmEN\n4HD/DnuamBzNfvf33NPM9Dh+efFsvqho4h/rDnYHO4DTwPdWbqK1z/7hyEmwuEgr09J8s6hKZkIU\n//WNudxz/kw+2F3NuzureG9nJS9uKu21XUxEGIkxEcTarESGhxEVHka41YLVIlgtQphFsIggQvcn\nRaH3p6W+TQ1/faIc7c6ckcpFcyf49RgjOlpGRG4CbgLIyckZ0j7q7Z1BtXbkQH/YItJ9XwQs4vqP\nEWYRrGGur21WC+OiwomItREXaSU+0kpcZDip8TbS4yPJGBfFxJQYYm2+/zXVuNfgdJhJAOwoF+4Y\nN4uUuCP99gnGYMvu4rODdWDgjH82YO9wsCBnGuFZTgpLGoixWbF3OliYm8iPvzGHxb97v/v10bkL\nicyxsrGojqhwC83t7jAzsPxdB8/eNA+LRXA0trHNHDnxZwG+lJfE7759ErX2TpJjI/oNiZqmdl6u\nOITDpAAwPzMBR9pcnj9Ug8M5nsJKV2Ctvus0blm5mc0l9QB8UQWzMiaz7dDRb1A936WjI8LA6cTe\ndeTBbZVHnrdZoIModpvWXruwAI9e9VUW/Olj9/ec6r710HH0ocH1N/TMt0/ipElJQ+5zH0hiTASX\nzM/kkvmZGGOobelgb1UzBw63cLi5ndqWTmpb2rF3OGjrctLW4aC1tZMup5Muh3E3UFwNA3D9qAZq\n4ADHdfJ5rJmeHuf3Y3iTGmVAdo/7We7HjncbjDErgBUA+fn5Q/rVnzEjlTNmpA6+oTqm5NgIFuYm\nuoIbWJibSHJs708EIsKfrlzA4gfew+E0NLW5WtWbi+v55MdnYrEISdHh1No7SYoO58rH1vV6/Q+e\nK+SfN55EXavr+V2VTVz4p49xGthUXE9NSwfj3a33L+UlsqmojgU5iTx89QJS4myISK83m/6+h/y8\nJDYerGVudgIvLD8ZEWFWRhxbyxoBaOlw0NjaxQvLT+Ybf/mUwpJ68vOS+MNlJ/KV33/Qa3+xtjCc\nxmDvcLWy27uc3UHWn3Yn7K5uPepxJ3D2gx/S2tn7tdHhFkSElo6jW+weX8pLZNHk8X5v8YoI42Nt\njI+1cZJe8BSSvAn3DcBUEZmIK7CXAlf12eZV4FZ3f/xJQIP2twc3EeHZm06murkdge4w7Ss13kZ+\nbiIbD9YSbbNi73C10lPjj2yfEmejuqmdTcX1vV67qaiOutbO7oCemRHPl/KSKCiqY2FuYvcbQkFR\nHQtyElj7k7N67deb7+GZ7yxynUzt0bp/+XuLmfvL1bS0O7q7NywWCy8sX9y97c7yxqP219Lu4OM7\nT+fsBz/C3uEg0mphRlosBSWuLqPoCEt38PdnzoRYtpU34zQcFewAbV3OXq3Z2ZnxPLZsYfenPM+b\nmXZlKF8YNNyNMV0icivwFq6hkE8aY7aLyHL3848Cq3ANg9yLayjk9f4rWfmKxSKkxUcec5ueAepp\npffXTeL5JNDzTWBBTqL7nINrKGPfMD7c7Ooa6nIaNhXXY7HIcQebxXJ0695qDaPwZ+ewt7q515BC\ni0UYH+M67rS0WKLCpVcIz86MJzw8jDZ3y7qlw4ET4dO7z6S+pYMLHvq413FibVbsHV3MzUrg0Wtc\nnzaufGw9G4vqiI4Io7mti1kZcRTV2mludxBjszIjLZZNxfXdnzR83fWilIcEav3F/Px8s3HjxoAc\nW/mHZ0hiUnQ4h5s7+P4zm9hUXM/C3ESe+c4igF6tbGMMS1es627JP3vTIr+2Wp1O0/1JYWFuIk9d\nm8+CX7+LvcNBjC2MLfeejcVi4bJH1nZ/CgmzCOvuPovk2AiWrljHxoO1zMkax6PXLCQ51nbUm13P\nn0GtvROHw8niB97DaSBMYO1PzsJikQHPIyg1GBEpMMbkD7adTj+gfKZnK9piETYV19PlNBQU1VHd\n1M4Pnt3cHayecev9dav4i+cksqempk4nW+49uoX//HdP5vK/fMqW0gby3eciBuoC6vupoefPYHxM\nBEtXfNo9+mlhbtJxdTspNRwa7sovPN00njAXoVew1rR0kBJn67dbZaRq8oR037HkYWEWXuzRP+8J\n4+Ottaal48gnAIE/XzVfg12NGA135Rd9W7rAUcEa6JqOFbS+eNPp+2YyUm9iSoH2uasR1D1NwBjq\nbx6L37PyL+1zV0FnJLtggsVY/J5VcNBxWEopFYI03JVSKgRpuCulVAjScFdKqRCk4a6UUiFIw10p\npUKQhrtSSoWggF3EJCLVQHAsqeSSDATzUvHBXh8Ef43BXh9ojb4Q7PXB8GrMNca9Qs0xBCzcg42I\nbPTmqq9ACfb6IPhrDPb6QGv0hWCvD0amRu2WUUqpEKThrpRSIUjD/YgVgS5gEMFeHwR/jcFeH2iN\nvhDs9cEI1Kh97kopFYK05a6UUiFIw70PEbldRIyIJAe6lr5E5L9EZKeIbBGRl0UkIdA1AYjIuSKy\nS0T2ishPAl1PXyKSLSJrROQLEdkuIj8MdE39EZEwEdksIq8Fupb+iEiCiLzg/hvcISInB7qmvkTk\n/7h/x9tE5BkROfYK8CNT05MiUiUi23o8liQiq0Vkj/vfRF8fV8O9BxHJBs4BigNdywBWA7ONMXOA\n3cDdAa4HEQkDHgbOA2YBV4rIrMBWdZQu4HZjzCxgEXBLENYI8ENgR6CLOIYHgTeNMTOAuQRZrSKS\nCfwAyDfGzAbCgKWBrQqAvwHn9nnsJ8C7xpipwLvu+z6l4d7bH4G7gKA8EWGMedsY0+W+uw7ICmQ9\nbl8G9hpj9htjOoBngSUBrqkXY0y5MWaT++smXKGUGdiqehORLOAC4PFA19IfERkHnAY8AWCM6TDG\n1Ae2qn5ZgSgRsQLRwKEA14Mx5kOgts/DS4Cn3F8/BVzi6+NquLuJyBKgzBhTGOhavHQD8Eagi8AV\nkiU97pcSZMHZk4jkAfOB9YGt5Cj/g6th4Qx0IQOYCFQDf3V3HT0uIjGBLqonY0wZ8Htcn7zLgQZj\nzNuBrWpAacaYcvfXFUCarw8wpsJdRN5x98X1vS0B7gHuDfIaPdv8FFdXw8rAVTr6iEgs8CLwI2NM\nY6Dr8RCRC4EqY0xBoGs5BiuwAHjEGDMfaMEPXQnD4e63XoLrjWgCECMi1wS2qsEZ15BFn/cWjKk1\nVI0xX+vvcRE5EdcfRKF7EeMsYJOIfNkYUzGCJQ5Yo4eIXAdcCJxlgmMcaxmQ3eN+lvuxoCIi4biC\nfaUx5qVA19PHKcDFInI+EAnEi8jTxphgCqZSoNQY4/nE8wJBFu7A14ADxphqABF5CVgMPB3QqvpX\nKSIZxphyEckAqnx9gDHVch+IMWarMSbVGJNnjMnD9Ye8YKSDfTAici6uj+4XG2Psga7HbQMwVUQm\nikgErhNYrwa4pl7E9Y79BLDDGPPfga6nL2PM3caYLPff3lLgvSALdtz/F0pEZLr7obOALwJYUn+K\ngUUiEu3+nZ9FkJ307eFV4Fr319cCr/j6AGOq5R4C/gzYgNXuTxjrjDHLA1mQMaZLRG4F3sI1OuFJ\nY8z2QNbUj1OAZcBWEfnc/dg9xphVAaxpNPo+sNL9Jr4fuD7A9fRijFkvIi8Am3B1W24mCK5WFZFn\ngNOBZBEpBe4Dfgv8S0RuxDU77jd9ftzg+GSvlFLKl7RbRimlQpCGu1JKhSANd6WUCkEa7kopFYI0\n3JVSKgRpuCulVAjScFdKqRCk4a6UUiHo/wOm0MSOmz//AQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f3858233ba8>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "samples = rejection_sampling(func, rv, k=15, n=100)\n",
    "plt.plot(x, func(x), label=r\"$\\tilde{p}(z)$\")\n",
    "plt.hist(samples, normed=True, alpha=0.2)\n",
    "plt.scatter(samples, np.random.normal(scale=.03, size=(100, 1)), s=5, label=\"samples\")\n",
    "plt.legend(fontsize=15)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 11.1.5重要性采样"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "回到期望的计算，我们考虑:\n",
    "$$\\begin{align}\\mathbb E[f]&=\\int f(z)p(z)dz\\\\&=\\int f(z)\\frac{p(z)}{q(z)}q(z)dz\\\\&=\\frac{1}{L}\\sum_{l=1}^L\\frac{p(z^{(l)})}{q(z^{(l)})}f(z^{(l)})\\end{align}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$r_l=\\frac{p(z^{(l)})}{q(z^{(l)})}$为重要性权重，修正了由于从错误的概率分布中采样引⼊\n",
    "的偏差。注意，与拒绝采样不同，所有⽣成的样本都被保留。\n",
    "* L个样本$z^{(l)},...,z^{(L)}$从$q(z)$中抽取\n",
    "* 计算权值 $$w_l=\\frac{\\tilde r_l}{\\sum_m\\tilde r_m}$$\n",
    "\n",
    "* L个样本的第二个集合从离散概率分布$z^{(l)},...,z^{(L)}$中抽取，概率由权值$(w_1,...,w_L)$给定"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAD8CAYAAACMwORRAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xd8XOWV8PHfmdGo17FkyVaxbNwBV0GMKQGzZOmGmICB\nOLSN16FssqEkIe+SDXkTkmx2s7DxS2JKIAsxYFpIKAEMprvbcsG9qPfey8zz/jEzsiRL1kiaJul8\nP5+xNXfu3Hskjc48c+65zxVjDEoppUYXS7ADUEop5Xua3JVSahTS5K6UUqOQJnellBqFNLkrpdQo\npMldKaVGIU3uSik1CmlyV0qpUUiTu1JKjUJhwdpxcnKyyc7ODtbulVJqRNq2bVulMSZloPWCltyz\ns7PZunVrsHavlFIjkojkebOelmWUUmoU0uSulFKjkCZ3pZQahTS5K6XUKDRgcheRSBHZLCK5IrJX\nRH7axzoiIo+JyGER2SUiC/wTrlJKKW940y3TBiwxxjSKiA34VETeNsZs7LbOZcA09+0rwOPu/5VS\nSgXBgMnduC7V1Oi+a3Pfel++aSnwJ/e6G0UkUUQmGGNKfBqtUqpP9fX1lJeX09HREexQ1DDYbDbG\njx9PfHz8sLflVZ+7iFiBbcBUYLUxZlOvVdKBgm73C93LNLmrUaukroXPDlfx9fnpWCwStDjq6+sp\nKysjPT2dqKgoRIIXixo6YwwtLS0UFRUBDDvBe5XcjTEOYJ6IJAKvicgZxpg9g92ZiKwEVgJkZWUN\n9ulKhYxdhbXc8exWKhraiLJZuWLOhKDFUl5eTnp6OtHR0UGLQQ2fiBAdHU16ejrFxcXDTu6D6pYx\nxtQCHwKX9nqoCMjsdj/Dvaz389cYY3KMMTkpKQOePatUSHpnTwnX/+ELwq0WJo2L5r/fP4jDGbwL\nzXd0dBAVFRW0/SvfioqK8kl5zZtumRT3iB0RiQIuAfb3Wu0N4FvurplFQJ3W29Vo9N6XZXzn+e3M\nmhDP63edy71fm8Gh8kbe3B3cl7uWYkYPX/0uvRm5TwA+FJFdwBbgPWPM30RklYiscq/zFnAUOAw8\nAdzpk+iUCjGvbCskLT6Std9eREpcBFecOYHpqbE8GuTRu1K9edMtswuY38fy33f72gB3+TY0pUJL\nh8PJZ4cruXLuBCJtVgCsFuF7/zCdO5/fzl9zi7lmfnqQo1TKRc9QVcpLO/JraWjr5KvTex4vuvT0\nNGamxfHY+kN0OpxBik6pnjS5K+Wljw6WY7UIi6cm91husQh3XjSVo5VNbMurCVJ0o8ctt9xCUlIS\njz76KGvXruXWW2/t8bgxhnnz5vHss896vc27776bO+64w8eRhjZN7kp56aODFSzISiQ+0nbSY+ee\nNg6AHQW1gQ5rVKmpqeHvf/87Tz75JM899xx33XUXt99+e491XnrpJaqrq7npppu83u59993H888/\nz+HDh30dcsjS5K6UFyoa2thTVH9SScZjXGwE2eOi2a4j92FJSkqitLSUZcuWsWXLFqqrq7ngggt6\nrPPYY4+xYsUKbLaT32T7k52dzXnnncfjjz/u65BDliZ3pbzw6eEKAL46fXy/68zPSmJHQS2u/gI1\nFHa7nX/7t3/j5ptvZty4cdjtdv75n/+ZtrY2AA4fPsznn3/Odddd1+N5L774IiLS562+vh6AZcuW\n8fzzz+N0jo3jIprclfLCRwcqGBcTzukT+z9rcEFWIhUNbRTVtgQwstEjPz+fmpoafvvb32K321m3\nbh3f+973WLNmDatXrwZg/fr1xMTEMHfu3B7PveCCC/jiiy+6bm+++SZ2u53LL7+860zPxYsXU1ZW\nxu7duwP+vQVD0K6hqtRI4XQaPj5UyVenp5xyDpn5WUkAbM+vJSNJpwIYrNzcXMB18POXv/wlAEuW\nLOHjjz/mnXfe4fvf/z7btm1j1qxZWCw9x6UTJkxgwgTXFBD19fVcfPHFnHHGGaxbt65rndNPPx2r\n1crmzZtPenMYjTS5KzWAPcV1VDe191tv95iRFkekzcKO/BqunjsxQNH17ad/3cuXxfVB2ffsifH8\n5KrTB/283NxcoqOjuf/++3ssnz59Ops3bwagtLSU5OTkvp4OQFNTE5dffjkWi4W//e1vPebbCQsL\nIzExkdLS0kHHNhJpcldqAJ8cqgTgvGn9JxUAm9XCnIxEduRrx8xQ5ObmcsEFFzBu3Lgey0tLS5k4\n0fVm2dra2u8Eaa2trVx99dU0NTXxwQcfEBcXd9I6ERERtLa2+j74EKTJXakB7CmqY3JyDMmxEQOu\nOz8rkac/PUZrh6PrLNZgGMrIOdhyc3M5//zzeyxraGjggw8+4Kc/dV0Azm639zny7ujoYNmyZZSU\nlPDRRx+RlJTU5z5qa2ux2+2+Dz4E6QFVpQZwoKyB6amxXq07PzOJDodhb5BKIiNVU1MTR44c4ciR\nIz2WP/bYY4DrxCaAGTNmcOzYsR7rOBwObrzxRg4cOMD7779PfzPOVlRU0NzczPTp0/3wHYQeHbkr\ndQqtHQ6OVzZx5Znezde+ICsRgB35NSyc1PfoUZ1s9+7dOJ1Ojhw5wsMPP8z555/Pu+++y3/+53/y\n3HPPkZjo+rmee+65PPzww1RUVHQl8e985zu89dZbPPPMM+Tn55Ofnw9ATEwMZ555Ztc+tm7dioiw\nePHiwH+DQaDJXalTOFLRiNPA9LST67d9GR8fSXpilNbdByk3N5f4+HjefPNNbrvtNn7+858zY8YM\n1q1bx9KlS7vWu/DCC7Hb7bzzzjusWLECYwwvvPACLS0t3HDDDT22edVVV/HGG2903X/nnXf46le/\nelJNf7TS5K7UKRwobQBgRqp3yR1cdXc9U3VwcnNzOfPMM5kzZw7btm3rd73w8HC++c1v8sILL7Bi\nxYoeJymdisPh4JVXXulqsRwLtOau1CkcKGvAZhWyk2O8fs6CrCSK61oprRsbXRm+kJub63Xv+f33\n38+HH37IwYMHvd7+unXriIqKYvny5UMNccTR5K7UKRwsbeC0lFhsVu//VOZmJgCuLhs1MGMMu3fv\nZs6cOV6tn5GRwdNPP01JifdXvzLG8NRTTxEWNnaKFWPnO1VqCA6WNZKTPbgDo1NTXCWco5WNQKof\nohpdvC2tdDfYEfiNN944qPVHAx25K9WPhtYOimpbmD6IejtAQrSNcTHhHK1o8lNkSg1Mk7tS/ThY\n1ggM7mCqx5SUGI5UNPo6JKW8psldqX50dcp42QbZ3WkpsTpyV0GlyV2pfhwsayA63Ep6YtSgnzsl\nJYaqpnZqm9v9EJlSA9PkrlQ/DpQ2MD017pTT/PZnSrJruoIjOnpXQaLJXal+HCxrGFK9HVwjd4Cj\nWndXQTJgcheRTBH5UES+FJG9IvLdPta5UETqRGSn+/aQf8JVKjAqG9uoamr3etqB3jLt0diswtFK\nHbmr4PCmz70TuNcYs11E4oBtIvKeMebLXut9Yoy50vchKhV4B4cw7UB3NquFLHu0jtxV0Aw4cjfG\nlBhjtru/bgD2Aen+DkypYNrvTu7T07yb6rcvU1JiteaugmZQNXcRyQbmA5v6eHixiOwSkbdFZORd\nKUCpbg6VN5AYbSPFiwt09Oe0lFjyqprodDh9GJlS3vE6uYtILPAK8D1jTO9zhbcDWcaYOcD/AK/3\ns42VIrJVRLZWVFQMNWal/O54ZTOTk2MQGXynjMeUlBg6HIbCmhYfRqaCZcOGDYgIe/bsCXYoXvEq\nuYuIDVdif94Y82rvx40x9caYRvfXbwE2ETnpgpPGmDXGmBxjTE5/V0tRKhTkVzczyd73tTq9dZqn\nY6ZS6+4q8LzplhHgKWCfMea/+lknzb0eInK2e7tVvgxUqUBp63RQXNdC1jjvp/ntS1eve7nW3VXg\neTNyPxdYASzp1up4uYisEpFV7nWuA/aISC7wGLDcGGP8FLNSflVY04IxDHvknhQTjj0mXEfuXtq7\ndy+XXnopdrudmJgYZs2axerVqwF48803ueSSSxg/fjzx8fEsWrSId999t8fz//3f/53k5GQ2bdpE\nTk4OUVFRnHfeeRw7dozy8nKuueYaYmNjmTVrFh988EGP52ZnZ3Pffffxs5/9jLS0NGJjY7n55pup\nqzv1tM1Op5Nf/vKXTJ06lYiICKZPn86zzz7bY51PP/2U888/n/j4eOLj45k3bx7r1q3zwU/s1AZs\nhTTGfAqcsvBojPkd8DtfBaVUMOVXNQOQnTy85A4wJTlGO2a8dNVVVzFr1iyee+45IiIiOHDgQNdU\nwMeOHePyyy/n3nvvxWq18vbbb3PZZZfx8ccfc+6553Zto7m5mZUrV/LAAw8QExPDv/zLv7BixQoi\nIiK47LLLuPPOO/n1r3/NN77xDQoKCoiOPvE7Xrt2LVOnTuWJJ56gpKSEBx54gH/6p386ZSK+5557\nePbZZ3nooYdYsGAB7733Hrfffjvjxo3jyiuvpL6+niuvvJKlS5fy0EMPdc1dX1sbgMswGmOCclu4\ncKFRKhT98dOjZtIP/mbK61uHva371+00C3/2ng+i6t+XX37p0+05HE5TXt9qnE6nT7d7KhUVFQYw\nu3btGnBdh8NhOjo6zNe+9jVz2223dS3/yU9+YgCzYcOGrmWrV682gPnpT3/atWzv3r0GMG+99VbX\nskmTJpmkpCTT0NDQtey5554zItL18/3www8NYHbv3m2MMebQoUNGRMwzzzzTI74VK1aYnJwcY4wx\nW7ZsMYCpr68fzI/jlL9TYKvxIsfq9ANK9ZJX3Ux0uJXk2PBhb2tKSiyVjW3UtXT4IDL/czoNNz6x\nkXMeWc/yNRtxOgNTXbXb7WRmZrJq1SpefPFFysvLezxeWFjILbfcQnp6OmFhYdhsNt59992TLrUX\nHh7O+eef33V/6tSpACxZsuSkZUVFRT2ee8kllxAbe+K8hmuvvRZjDFu2bOkz5vXr12OxWLj22mvp\n7Ozsul188cXs3LkTh8PBaaedRmxsLDfddBN/+ctfAjNid9PkrlQv+VXNZNmjh9UG6XFaiitZjJQz\nVaua2tmWV0On07Atr4aqpsDMammxWHj33XdJS0vj9ttvJy0tjfPPP58dO3bgdDq5+uqr+fzzz3n4\n4Yf58MMP2bJlC5dddhmtrT2vUxsXF4fFciKthYe73qATExNPWtb7uePHj+9xPzo6mtjY2H4v51dZ\nWYnD4SAhIQGbzdZ1u/XWW+ns7KSkpISkpCTee+89Ojo6uP7660lJSeGKK67g6NGjQ/9heUkvs6dU\nL3nVzV1tjMM12X1h7WOVTczPGtzl+oIhOTachZOS2JZXw8JJST759OKtmTNn8sorr9DR0cEnn3zC\nD37wA6644go2bNjAjh07ePvtt7n00ku71m9p8e35A70/LTQ3N9PY2MiECRP6XN9utxMWFsZnn33W\n4w3Fw/NmsWjRIt555x1aWlp4//33+f73v89NN93Exo0bfRp/b5rclerG6TTkVzezZOb4gVf2QkaS\nay74guqRcSKTiLD224uoamonOTbcJ59eBstms7FkyZKuJOgZOUdEnDhbOC8vj88++8zri2p74733\n3qOxsbGrNPPaa68hIuTk5PS5/pIlS3A4HNTV1XHJJZcMuP2oqCiuuuoq9uzZwyOPPOKzuPujyV2p\nbsoaWmnvdJI1zDZIj0ibldT4CApqmn2yvUCwWISUuKFPuzAUu3bt4r777uOGG25gypQp1NTU8Ktf\n/Yq5c+eyaNEiMjIyuPfee/nZz35GQ0MDP/nJT0hP9+0UV1FRUVxxxRXcf//9lJSUcP/993Pttdcy\ne/bsPtefMWMGq1atYvny5TzwwAPk5OTQ2trK3r17OXjwIE8++SRvvvkmTz/9NNdccw1ZWVkUFRXx\nhz/8occxAH/R5K5UN3nuNshJ43yT3AEyk6IpqB45yT0Y0tLSSE1N5ec//znFxcUkJiZy0UUX8atf\n/YqIiAheffVV7rrrLq677joyMjL48Y9/zIYNG3w6FcDy5cuJi4vjjjvuoLGxkauvvprHH3/8lM9Z\nvXo106dP54knnuChhx4iPj6e2bNnc8cddwCug7ciwoMPPkh5eTkpKSlceeWV/OIXv/BZ3P0RE6Rz\njXJycszWrVuDsm+l+vPSlgIeeGUXH99/EVk+SvD/+uJONh+r5rMf+me0tm/fPmbNmuWXbY8V2dnZ\nXHfddfzmN78JdijAqX+nIrLNGNN3ragb7ZZRqpu86ibCLMLExEifbTMzKYqSuhbaO3V2SBU4mtyV\n6iavqpn0pCjCrL7708i0R+M0UFw7Mg6qqtFBa+5KdZNf3eyzg6keme7tFdQ0k53smxZL5VvHjx8P\ndgg+pyN3pbrJq2r26cFU6JbcR0g7pBodNLkr5VbX3EFdSweT7L4dXafFR2Kzil/bIYPVGKF8z1e/\nS03uSrnlVbtmb/RVl4yH1SJMTIzyWzukzWbz+dmaKnhaWlqw2WzD3o4md6Xc/NHj7pFlj6bAT5fb\nGz9+PEVFRTQ3N+sIfgQzxtDc3ExRUdFJ89wMhR5QVcot3z2y9vUBVYCMpGj+vrfU59sFiI+PB6C4\nuJiOjpEx+6Tqm81mIzU1tet3Ohya3JVyy6tqYnxcBNHhvv+zyLRHUd3UTlNbJzERvt++5yo/Snlo\nWUYpt4Lqlq7OFl/LTDrRDqlUIGhyV8qtoKaZTPcsjr6m7ZAq0DS5KwV0OpyU1LWSkeSvkbtn6l8d\nuavA0OSuFFBa34rDabrmX/c1e0w4MeHWroO2SvmbJnelgEJ3m6K/Ru4iQqY9mkKtuasA0eSuFN2T\nu39G7q5tR2vNXQXMgMldRDJF5EMR+VJE9orId/tYR0TkMRE5LCK7RGSBf8JVyj8Ka5oRgQk+nOq3\nt0x7FAU1eqKRCgxvRu6dwL3GmNnAIuAuEel93anLgGnu20rg1JcvUSrEFNa0kBYfSUSY1W/7yEyK\nprndQXVTu9/2oZTHgMndGFNijNnu/roB2Af0vnjhUuBPxmUjkCgifV8yXKkQVFjT7NeSDJw481UP\nqqpAGFTNXUSygfnApl4PpQMF3e4XcvIbgFIhq7CmxW8HUz1OzOuudXflf14ndxGJBV4BvmeMqR/K\nzkRkpYhsFZGtFRUVQ9mEUj53osfdvyN3z/a1Y0YFglfJXURsuBL788aYV/tYpQjI7HY/w72sB2PM\nGmNMjjEmJyUlZSjxKuVz/u5x94iJCCMp2tbVmaOUP3nTLSPAU8A+Y8x/9bPaG8C33F0zi4A6Y0yJ\nD+NUym887Yn+Lst49lGkyV0FgDfT050LrAB2i8hO97IHgSwAY8zvgbeAy4HDQDNwm+9DVco/PGUS\nf4/cAdITozhU3uD3/Sg1YHI3xnwKyADrGOAuXwWlVCAV1rS4etwT/J/cM5Ki2HCwHGMMrg/FSvmH\nnqGqxjxPj3t4mP//HDKSomjtcFKlve7KzzS5qzEvED3uHunuur4eVFX+psldjXmB6HH38LyJ6EFV\n5W+a3NWY1ulwUlrv/x53j3TtdVcBosldjWkldYHpcfeIj7QRHxmmZRnld5rc1Zjm73nc+5KRFE1R\nrSZ35V+a3NWYFsged4+MpCgtyyi/0+SuxrRA9rh7pCdFUVjTovO6K7/S5K7GtED2uHtkuOd1r23u\nCNg+1dijyV2NaQUB7HH3ODE7pNbdlf9ocldjWlFNC5kBPJgKrvllQNshlX9pcldjVofDSUldS8BH\n7p43E+2YUf6kyV2NWcW1LTgNZNgDO3KPjwojLkJ73ZV/aXJXY5YnuQa6LCMiXR0zSvmLJnc1ZhVU\nB77H3UN73ZW/aXJXY1ZBTTNWizAhITLg+05PjNLJw5RfaXJXY1ZhTQsTEyMJswb+zyAjKZqGtk7q\nWrTXXfmHJnc1ZhVUN5ORGNh6u0eGzg6p/EyTuxqzCmpayLQHvt4O3af+1dKM8g9N7mpMau1wUNHQ\nFvBOGY8MvSKT8jNN7mpM6prqN0gj96RoG9HhVi3LKL/R5K7GpAJ3Ug3WyF1EyEyKpqBaR+7KPzS5\nqzGp6wSmAJ+d2l2mXXvdlf8MmNxF5GkRKReRPf08fqGI1InITvftId+HqZRvFVY3Ex5mISU2Imgx\nZCRFk1/drPO6K7/wZuT+DHDpAOt8YoyZ5749PPywlPKvgppmMhKjsFgkaDFk2l3zulc3tQctBjV6\nDZjcjTEfA9UBiEWpgCmsaQn4hGG9ZbrbIQu0Y0b5ga9q7otFZJeIvC0ip/tom0r5TUF1c1dyDZas\ncdFdsSjla2E+2MZ2IMsY0ygilwOvA9P6WlFEVgIrAbKysnywa6UGr7Gtk5rmjq5e82DxdOoU6EFV\n5QfDHrkbY+qNMY3ur98CbCKS3M+6a4wxOcaYnJSUlOHuWqkh8XSoBOvsVI+YiDDsMeHaDqn8YtjJ\nXUTSRETcX5/t3mbVcLerlL94kmmwety7y9Spf5WfDFiWEZG1wIVAsogUAj8BbADGmN8D1wHfEZFO\noAVYbrS3S4WwYM7j3luGPZq9RXXBDkONQgMmd2PMjQM8/jvgdz6LSCk/K6xpITrcij0mPNihkJkU\nzbt7S3E4DdYgtmWq0UfPUFVjTkFNM5lJ0biriUGVaY+iw2Eoq28NdihqlNHkrsacgurmkCjJwIm6\nf762Qyof0+SuxhRjDPnVzUwaFxPsUIATc9tor7vyNU3uakypaGijud1BdnLwO2UAJiZGIqJnqSrf\n0+SuxpTjVa4RcqiM3CPCrKTFR1KoI3flY5rc1ZhyvKoJgOxxoTFyB1fdXc9SVb6myV2NKccrmwiz\nCOmJoXFAFVx1dz1LVfmaJnc1puRVNZNpjybMGjov/Ux7FGUNrbR1OoIdihpFQucVrlQAHK9qYlII\nlWTAVZYxBor0oKryIU3uaswwxnC8sonsEDmY6tHVDqnJXfmQJnc1ZlQ2ttPU7gipg6lwYnZK7XVX\nvqTJXY0Zee5OmUnJoTVyT42LJNxq0eSufEqTuxozjlV62iBDK7lbLEKmPaqrTVMpX9DkrsaMvKpm\nrBYJmXllupucHNv15qOUL2hyV2PG8aomMpKisIVQG6TH5ORojlc143TqpRCUb4Teq1wpP8mrCp0J\nw3rLTo6hvdNJcZ12zCjf0OSuxoQTbZCh1SnjMdl9kPd4pR5UVb6hyV2NCdVN7TS0dYbcwVQPT3I/\npgdVlY9ocldjgmc2yFCZ6re31LhIomxWjlVocle+ocldjQnH3Z0ooVpzt1iESeOitR1S+YwmdzUm\n5FU1YZETl7ULRZOTY7rehJQaLk3uakw4XtVMelIU4WGh+5KfnBxDfnUznQ5nsENRo0DovtKV8qHj\nVaE3YVhv2ckxdDoNhTqBmPKBAZO7iDwtIuUisqefx0VEHhORwyKyS0QW+D5MpYbOGMOxitCb6re3\nKdoxo3zIm5H7M8Clp3j8MmCa+7YSeHz4YSnlO6X1rTS0dTIjNS7YoZxStie5a8eM8oEBk7sx5mOg\n+hSrLAX+ZFw2AokiMsFXASo1XAdKGwCYHuLJfVxMOHERYdoxo3zCFzX3dKCg2/1C9zI1RlQ2trH6\nw8P8eVM+rR2hd6m4g2UjI7mLCJNTYnQCMeUTYYHcmYisxFW6ISsrK5C7Vn5QXNvCmo+PsnZzPm2d\nrg6P/3z3ALcuzuZbi7NJiLIFOUKXA6WNjI+LICkmPNihDCh7XAzb82uCHYYaBXwxci8CMrvdz3Av\nO4kxZo0xJscYk5OSkuKDXatgKahu5rJHP+G5jXlcPXci6+/9Ki+sXMScjAT+872DLF+zMWRG8QfL\nGpiRFtqjdo/JyTEU1bboxbLVsPkiub8BfMvdNbMIqDPGlPhguypEdTic3LN2B05jeOd7F/Af35jL\naSmxLJoyjj/edjZPfiuHfSX1/PSvXwY7VBxOw6HyhpAvyXhMTo7BGMiv0gnE1PB40wq5FvgCmCEi\nhSJyh4isEpFV7lXeAo4Ch4EngDv9Fq0KCb959wA7C2r51bI5TB0fe9Lj/zA7le9ceBprN+fzl519\nfogLmILqZlo7nCHfKePRNYGY1t3VMA1YczfG3DjA4wa4y2cRqZC24UA5f/joKDd/JYvLz+y/Kere\nS6az5Vg1D766mzPSEzgt5eQ3gUA44DmYOkLKMp52SO2YUcOlZ6gqrzW0dnDfulxmpsXxb1fOPuW6\nYVYL/3PTfCJsVu59KRfXGCDwDrrbIKf18QkjFCVE2bDHhHNUe93VMGlyV157bmM+lY3t/GrZHCJt\n1gHXn5AQxQP/OIOdBbVsOFgRgAhPdqCsgUx7FDERAW0MG5Zp42O72jeVGipN7sorze2dPPnJUb46\nPYW5mYleP+/rCzJIT4ziv98/FJTR+8GyhhFTb/eYmRbHwbJGvZ6qGhZN7sorazcXUNXUzj1Lpg7q\neeFhFu5eMpXcglo+CvDovb3TydGKphHTKeMxIy2exrZOimp1AjE1dJrc1YBaOxys+fgIi6bYycm2\nD/r5y9yj90fXB3b0fryqiU6nGTE97h6eePeXamlGDZ0mdzWgl7cVUlbfxj1Lpg3p+eFhFu66aCo7\n8mv5+FClj6Pr30iZU6Y3T3I/UFof5EjUSKbJXZ1Sh8PJ4xuOMD8rkcWnjRvydq5b6Bq9/8/6Qz6M\n7tQOljVgtQhTUkJ7HvfeYiPCyLRH6chdDYsmd3VK6/eVU1Tbwne+ehoiMuTthIdZuOO8yWzNq2FP\nUZ0PI+zfgdIGJifHEBE2cGdPqJmRGq/JXQ2LJnd1Sms355MWH8mSmeOHva1lCzOIsll5bmOeDyIb\n2EjslPGYmRbHscomnWNGDZkmd9WvgupmPj5UwfVnZRJmHf5LJSHKxjXzJ/L6ziLqWjp8EGH/6ls7\nyKtuZuYIO5jqMSMtDofTcLi8MdihqBFKk7vq10tbXdP033BW5gBreu+biybR2uHk5W2FPttmX3bk\n12IMLJiU5Nf9+MvMroOqWppRQ6PJXfWp0+Hkpa0FXDg9hfTEKJ9t9/SJCSzISuS5jXl+PUln2/Fq\nLALzBnHCVSiZnBxDuNWiyV0NmSZ31acPD1RQVt/GjWf7/qIqK86ZxLHKJj4/UuXzbXtsy69h1oT4\nETXtQHfCfyjqAAAYrUlEQVRhVgtTx8eyT5O7GiJN7qpPazfnMz4uwicHUnu7/MwJ2GPC+dMXx32+\nbXB96tiRX0vOCC3JeMxMi9NedzVkmtzVSUrqWthwoJzrc3xzILW3iDArN5yVyfr95ZTWtfp8+/tL\nG2hud4zYervHjLQ4yurbqG1uD3YoagTS5K5O8ur2IpwGrs/x3YHU3m48KwuH0/DiloKBVx6kbXmu\na5AOZaqEUKLTEKjh0OSuejDG8Mr2Qs7OtpM1Ltpv+8kaF83505J5cUs+Dh8fWN2aV8OEhEifHggO\nhlkT4gHYX6KlGTV4mtxVDzsLajla0cSyhel+39fNX8miuK6VDQfKfbrdbcerR3xJBmB8XASJ0bau\nq0kpNRia3FUPr2wvJNJmOeUl9Hzl4lmppMRF8OdN+T7bZnFtC8V1rSP+YCqAiDB7Qjy7CgMzXYMa\nXTS5qy5tnQ7+mlvCP56eRlykze/7s1ktXJ+TwYcHyin20dzlnnr7wlGQ3MH1fewvbaCprTPYoagR\nRpO76rJ+Xzl1LR0sW5ARsH0uPysLA7zgowOr2/JqiLJZu+rVI92CSUk4nIbcgtpgh6JGGE3uqssr\n2wpJi4/k3KnJAdtnpj2aC6al8OKWfDoczmFvb1teDfMyE7H5oYUzGBZkuj6BeD6RKOWt0fEXoIat\noqGNDQcruGZ+OlbL0Kf2HYpbFk+irL6Nt3aXDGs79a0dfFlSP2pKMgAJ0Tamp8ayLV+TuxocTe4K\ngNd3FOFwGq4LQJdMbxdOH8+UlBie+vTYsC7D9+7eMhxOw5JZvj+rNpgWTkpie16NXjBbDYpXyV1E\nLhWRAyJyWER+2MfjF4pInYjsdN8e8n2oyl+MMby0tYD5WYlMHR/4KXItFuGO8yazq7COzceqh7yd\nN3KLybRHMX+EThbWnwVZSdS3dnK4Qqf/Vd4bMLmLiBVYDVwGzAZuFJHZfaz6iTFmnvv2sI/jVH6U\nW1jHofJGv56ROpBlCzJIirbx5KfHhvT8ysY2PjtcyVVzJg7rilGhyFNm0rq7GgxvRu5nA4eNMUeN\nMe3AC8BS/4alAumlrQVE2ixcOcf/ve39ibRZWbFoEu/vK+NYZdOgn//W7hIcTsPV8yb6Ibrgmpwc\ngz0mXJO7GhRvkns60L1PrdC9rLfFIrJLRN4WkdP72pCIrBSRrSKytaKiYgjhKl9raXfw153FXH7m\nhID0tp/KN8+ZhM1i4ekhjN7f2FnMjNQ4ZqaNjhbI7kSEBVmuurtS3vLVAdXtQJYxZg7wP8Drfa1k\njFljjMkxxuSkpKT4aNdqON7ZW0JDW2dQSzIe4+MiWTpvIuu2FVBe7/1skUW1LWzNqxmVo3aPhZOS\nOFrZRHWTzhCpvONNci8Cuv/lZ7iXdTHG1BtjGt1fvwXYRCRwzdJqyF7aUkiWPZqvTA6NGRTvumgq\nTif88u39Xj/nr7nFAFw1Z3Qnd9C6u/KeN8l9CzBNRCaLSDiwHHij+woikibuo1gicrZ7u/67zI7y\niYLqZr44WsU3FmaEzEHI7OQYVl4whVd3FLHluHedM3/ZWcy8zES/zmIZbHMyErBZRZO78tqAyd0Y\n0wncDfwd2Ae8ZIzZKyKrRGSVe7XrgD0ikgs8Biw3w2lYVgHx/KZ8LALLFgZuugFv3HnRaUxMiOSh\nv+wdcDrgzw5Xsq+knqWjuCQDrgPOp09MYPMxHTMp73hVczfGvGWMmW6MOc0Y83P3st8bY37v/vp3\nxpjTjTFzjTGLjDGf+zNoNXytHQ5e2JLP12anMTHE5j2PDg/j/1w5m30l9fx5U16/69W1dHDfulym\npMSw/CzfX+s11Hx1ego7C2qpbGwLdihqBNAzVMeoN3YWU9vcwS2Ls4MdSp8uOyONxaeN4zfvHuRQ\nP/OZ//SNvZQ3tPHb6+cRFW4NcISBd8nsVJwGPtjn2/nv1eikyX0MMsbwzOfHmZEax6IpoXEgtTcR\n4f9ecwbhYRa+/vjnfH6kssfjb+8u4dUdRdx90VTmjrIzUvtz+sR40hOjePfLsmCHokaAsGAHoAJv\na14NX5bU84trzwyZA6l9mZISy2t3Lub2Z7Zwy9Ob+bcrZxMXGcbeonpe3l7InIwE7l4ydXg7Kd7h\n/boT5w9vX8MkIlwyO5W1m/Npbu8kOlz/fFX/dOQ+Bj37+XHiI8O4Zn7oH4TMSIpm3arFnJVt56G/\n7OVfX8zlfzfmMTUllv++Yd6omdrXW1+bnUpbp5OPD1YOvLIa0/Stf4wpq2/lnT2l3Lo4e8SM/BKi\nbDxz29lsPFpFWkIkU5JjCBtjSd3jrMl24iPDeO/LMi49Iy3Y4agQNjL+upXPPP3pMZzGsOKcScEO\nZVDCwyxcMF3ParZZLVw8K5X1+8vodDjH7JucGpgm9zGkoqGNZ784zjXz0pk0LibY4Yw8IVKfv2R2\nKq/tKGLL8RrOOW2c3/ajRjZ92x9Dfv/RETochnsunhbsUNQwXDA9hfAwC+9p14w6BU3uY0RZfSvP\nbczj6/PTmZyso/aRLDYijPOmJvPW7hKfXHdWjU5alhkjHt9wBIfTcM8SHbUHhJ9LODedncU//Wkr\nb+0uYem8wF8aUYU+HbmPASV1Lfx5Uz7XLcwY1ZNrjSVLZrquO/vkJ8O77qwavTS5jwG/fHs/BsNd\nFw3zhB8VMjzXnd1dVMemYVx3Vo1emtxHuQ/2l/GXncXcddFUMu06ah9Nli3IwB4TzpOfHA12KCoE\naXIfxRpaO/jxa3uYnhrLnRfqqH20ibRZ+eaiSby/r5wjFY3BDkeFGD2gOor9+p0DlNa3svrmxYSH\njaH38cEczAwFwzj4+q1zJvH7j47w1KfH+MW1Z/o4MDWSjaG/+LFl09Eq/ndjHrctnsyCrKRgh6P8\nJDk2gmULMli3tYB9JfXBDkeFEE3uo1BBdTN3/Xk7k8ZFc98/Tg92OMrP7vvadBKibHz/pVzaO7Xv\nXbloWWaUqW1u55Y/bqbDYXjqlrN8NzlYsE+9H2mllgAaFxvBI1+fw7f/tJVH1x/k/n+cGeyQVAjQ\nkfso0tbpYOX/bqOwuoU1KxYydXxssENSAXLJ7FS+sTCDxzcc0YtoK0BH7qNGXUsH331hB5uPVfPo\n8nl8ZYoXE0r5azQc7FF+L04n1LV2kBhlo69rkwz0+Ejx0FWz+fxIFf/64k5e/OdFTEgIrWvjqsDS\nkfsocKisgWtWf8anhyp55Otn6uno3Tid8ODru7n1j5v50Wu7cToH93iwOJ1Q09zBYE4+jYu08diN\n86luamfZ//ucw+V9X3tWjQ2a3Ecwp9Pw2o5Crln9GQ2tHfz524u48eysYId1kr4S1Yllvjt1vq/9\n1LV2sK+kHofTsK+knrrWjh7PGehxX+l0GI5VNuF0Dvz9DucNZ+GkJF5YuYh2h+G633+hJZoxTMsy\nI5Axho8OVvDrdw7wZUk98zITefybC4LyMdybkseDr+9mX3Ed01Lj+NWyOQjiWlZST2SWg7XfXoTF\nMrh6SO/9djoMP3h1F4fLG5k1IZ5fXHMmFgskRtmYNSGefSX1TB0fS0Jkz5d8YpSNmWnx7C+tZ+aE\neBKjbF5/n96Wc9raHdz89GZaOxzERITx/B1fIcza/xPqWjvYV1yHw8C+4jrqWjtIiu47rr6ckZ7A\nq99ZzIqnN3Hjmo1865xJ3HXRVJJiwr3ehhr5vBq5i8ilInJARA6LyA/7eFxE5DH347tEZIHvQ1WF\nNc2s+fgIS1d/xq1/3EJDWwf/fcM8Xv3O4iEn9qF8/O/+3B+9totbntrEfet2UtnQdtJ2uieq/aUN\nPPDKLmqa27tGy9vyaqhqah/0fruPbDs6nNz7ci4HShu6RuB51U0YAyLwf5eewbTxsRwua+DB1/f0\nGAm74jUYp6Gz04Gz1zfgdEJZXSv/8sI2vvXUpq6RdKfDcP8rudzy9CbufzmXjk5nnz/HTofh5qc3\n0drhAKCprZP86iaqmtqpamyjuqkdY3r+HuIjw4h0dzlFhocRHzn4MVjWuGhe+c5irp43kac/O8YF\nv/6QR98/RF5V06C3pUamAV81ImIFVgOXAIXAFhF5wxjzZbfVLgOmuW9fAR53/6+GoMPhpLSulZK6\nVvaX1rOnqI5dhXXsL3XVUM9MT+BnS0/nhrOyhnXmadeouqS+x2i39zp9jVhrmtupbW5nb7HrxJkD\nZY3c+swWTp8YzyPXzunaTmKUjWmpcV2xHy5vBEzXaHphVhLJseE4nYaqpnaSY8ORAY5q9i6l3P/q\nLo6Unzj9PjzMwvde2MGsiQn84pozaWjr5FB5Iw4DXxbVsaeoltkT42loc+BwOtlXXI/T/T088HIu\nD14+G3tMOMbAD1/J5cvSE7XrPUV1HKts5P9tOMyBMtc+95c2cPNTm2jrcDBzQjwP/ONM7DHhiEBB\nTTOtHSfeTSwCj284zL7SE/FOGx+DRYTDFU3u58/oejNo7XBQ39o5qJG7R3JsBL/5xlxWXjCF//j7\nAX77/kF++/5Bpo6P5aIZKcyeGM/UlDimpMQQE6Ef4kcbb36jZwOHjTFHAUTkBWAp0D25LwX+ZFwF\n1I0ikigiE4wxJb4OeG9xHa9uL/L1ZofNGDCYrq9PLHctdRqD07jq5J1OQ6fDSbvDSVuHk9ZOB42t\nnTS0dlLf2kGVezTnYY8J54z0BJbOS+eKMyf4bNreU3389yTwX/99P/tLG7qSP7hG656k3tt+d93a\nsx0R+NWyOTzw8i4OlLlG1r9+Zz8/v3YODW2dJEbZMEU7B3yT6S4+Moyp42M5XN7I5OSYHok9e1w0\n+dXNOA3sK6nneFUjCVE2Zk6IZ19RHSLw4Ot7sAo4DUSHh9G9pH2grJHbn93CzLQ4Vp4/pUdi9/ju\niztPWtbc7krGe4vruf3ZLUwdH8uvl81hkj2a6HBr1+PGwP7SnvPAHCo/MZreW1QH5sSb36xTlIq8\nNT01jie+lUN+VTPv7yvjg/3lPPP5cTocJ15kMeFW7LHh2KPDiYkII8pmJTLcSrjVgtUi2KyCRTw3\nut6ARUA48WY8kruNAmnxaeO4eFaqX/fhTXJPBwq63S/k5FF5X+ukAz2Su4isBFYCZGUN7cBfQXUL\nL24pGHjFIJCuf1z/ef4ALAIWEUTAahHCLK4/mIgwCxE2C+FWCwnR4WTao4mLtDE+LoIJCZFMSIxi\n6vhYJiZEDjiaHYr4yTkcDa+hobWTuIgw4ifngNVCZ6eTb6z5gl2FdTicyUAye4uFexNm43QaXiwq\nBpL73ObZk5JIPO2sHn/lVuC+22ZxziMfALCnGL6fMJvU+EgAKhvaeLm4hE7nuK79pMRF9Ll9p9Nw\n0xMb2VySQlR4KrllTmIjJtDU5mBeZgKPrjqHm57czNa8GqLDrCxZ63oTysmawj/fOIW7/uxu0/Tk\ntbY+duKA3CJ4+cUGHGbKIH+q7ueXwIHXmvjLnefy/P+Zy7W//5wvi+vJmZTE3uI6mtr7P0raMX4u\nv7jzK15/kvFW1rhobj9vMrefN5n2Tid5VU0cLm/kaGUTlY1t1DS1U9XUTku7g7qWDlo6HHQ6XAOR\nDqfBuAcoDvfX7opWF51V3nsxEWEhkdx9xhizBlgDkJOTM6TXwqVnpHHpGWk+jWusqm7u6BpRNrc7\nqG7uYFxMONev+YId+bVd61nF1YWRHBtORcPJ2fDNe84lOTYCESElLqLPZGTttaz7veTYcBZOSmJb\nXk3XfvpT1dTO1uPVGKDZnSAb2xxYRbCFWRGxsPbbizhY1sDlj33S9byt+XVs/XPP/nvh1AnJMcgD\nEQLMnhjH3mLXaH9PUT1X/+4T3rj7fN646zz3JzLDOY+s73rO6RNi2VvScyTf38HlwZSuBhIeZmFa\nahzTUuOGtR0Vurwp2BYBmd3uZ7iXDXYdFWKSY8PJmZREmEXIcSfVqqZ2cgtOJPb5mYl88aOLeWHl\noq7kfVb2iYnIzspOYvbEBFITohgff+IThtNpqGho62p1TImL4OxsO1aLcPZke4+RuYiw9tuLeuzn\nVDHPzUw8abnDGLa7D85aLMKMtDhysu19bkOAnKwENv5oCQuzEk75M4oJt/a4bxXhjIl9J0QD5FU1\n91i2t6SRZY9/Drh+BilxEeRk2wmzCF+ZbOeNu89n44+WMC8jHgHmZyZgj7Zx4xMbOeeR9SxfsxGn\n0+B0mpOWKXUqMlCfsYiEAQeBi3El7C3ATcaYvd3WuQK4G7gcV8nmMWPM2afabk5Ojtm6devwoldD\n5hkF2qNtVDd3dI0GjTEsX7ORrcermZuZyMurzsHSqwDuSdwi9DlS9yQiz0jc0+roq5Gnw+Hkuj+4\nykZnpscTbrWwPb+2q8e79xsMGO55YSfb82pYkJXI725a0BW302koq2+lurGNH7++h11FdXjyplXg\nsx9cxKrntrOzsA5wvdntKqzFYVzlNm9yrEVg04P/0PWG1tfPoXspbE5GArsKXPsIswhf/OhiAM55\nZD2dTtO1rL/SlRrdRGSbMSZnoPUGLMsYYzpF5G7g77jKp08bY/aKyCr3478H3sKV2A8DzcBtwwle\n+VdfyffEATLXKPpUSdhiEVITIvvdflVTO9vyaujs1uqYEheBxSI+SUhWq4VXVi3uitEY+oy3e5wv\n9PM9WSzChMQoJiRG8eqd51LZ2Mbda3ew3f2zSU1wLa9obENwfXK48YlNbMurYcGkJDCGbfm1RIdb\naWrtJCYyjOZ2R1d9Gjip1NTXz6GmpYPdhXU4nIZdhXXMzUxkV2Fdj+d6W7pSCrwYufuLjtyDp6Kh\nza+jQM/o35OIBiq1hJqBPmF0f9zzxuL5BOT5PzEyjEMVjYyLCe9RrupP75/Zn//pKz0+UXkTlxob\nvB25a3IfgwKRfDURDZ7+zJQ3fFaWUaOPN6WX4fJVCWYs0Z+Z8iVN7mOUJhKlRjedFVIppUYhTe5K\nKTUKaXJXSqlRSJO7UkqNQprclVJqFNLkrpRSo5Amd6WUGoWCdoaqiFQAeUHZed+SgcpgB3EKoR4f\nhH6MoR4faIy+EOrxwfBinGSMSRlopaAl91AjIlu9OaU3WEI9Pgj9GEM9PtAYfSHU44PAxKhlGaWU\nGoU0uSul1Cikyf2ENcEOYAChHh+EfoyhHh9ojL4Q6vFBAGLUmrtSSo1COnJXSqlRSJN7LyJyr4gY\nEUkOdiy9ich/iMh+EdklIq+JyMlXig4CEblURA6IyGER+WGw4+lNRDJF5EMR+VJE9orId4MdU19E\nxCoiO0Tkb8GOpS8ikigiL7tfg/tE5Jxgx9SbiPyr+3e8R0TWikj/14MMXExPi0i5iOzptswuIu+J\nyCH3/0mn2sZQaHLvRkQyga8B+cGOpR/vAWcYY+bgumj5j4IcDyJiBVYDlwGzgRtFZHZwozpJJ3Cv\nMWY2sAi4KwRjBPgusC/YQZzCo8A7xpiZwFxCLFYRSQf+BcgxxpyB65rPy4MbFQDPAJf2WvZDYL0x\nZhqw3n3fpzS59/Rb4AEgJA9EGGPeNcZ0uu9uBDKCGY/b2cBhY8xRY0w78AKwNMgx9WCMKTHGbHd/\n3YArKaUHN6qeRCQDuAJ4Mtix9EVEEoALgKcAjDHtxpja4EbVpzAgSkTCgGigOMjxYIz5GKjutXgp\n8Kz762eBa3y9X03ubiKyFCgyxuQGOxYv3Q68HewgcCXJgm73CwmxxNmdiGQD84FNwY3kJP+Na2Dh\nDHYg/ZgMVAB/dJeOnhSRmGAH1Z0xpgj4Da5P3iVAnTHm3eBG1a9UY0yJ++tSINXXOxhTyV1E3nfX\n4nrflgIPAg+FeIyedX6Mq9TwfPAiHXlEJBZ4BfieMaY+2PF4iMiVQLkxZluwYzmFMGAB8LgxZj7Q\nhB9KCcPhrlsvxfVGNBGIEZFvBjeqgRlXy6LPqwVj6hqqxph/6Gu5iJyJ6wWR675YdAawXUTONsaU\nBjDEfmP0EJFbgSuBi01o9LEWAZnd7me4l4UUEbHhSuzPG2NeDXY8vZwLXC0ilwORQLyIPGeMCaXE\nVAgUGmM8n3heJsSSO/APwDFjTAWAiLwKLAaeC2pUfSsTkQnGmBIRmQCU+3oHY2rk3h9jzG5jzHhj\nTLYxJhvXC3lBoBP7QETkUlwf3a82xjQHOx63LcA0EZksIuG4DmC9EeSYehDXO/ZTwD5jzH8FO57e\njDE/MsZkuF97y4EPQiyx4/5bKBCRGe5FFwNfBjGkvuQDi0Qk2v07v5gQO+jbzRvALe6vbwH+4usd\njKmR+yjwOyACeM/9CWOjMWZVMAMyxnSKyN3A33F1JzxtjNkbzJj6cC6wAtgtIjvdyx40xrwVxJhG\nonuA591v4keB24IcTw/GmE0i8jKwHVfZcgchcLaqiKwFLgSSRaQQ+AnwS+AlEbkD1+y41/t8v6Hx\nyV4ppZQvaVlGKaVGIU3uSik1CmlyV0qpUUiTu1JKjUKa3JVSahTS5K6UUqOQJnellBqFNLkrpdQo\n9P8BOEpKRKvzoo0AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f34e8ea32e8>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "samples = sir(func, rv, n=100)\n",
    "plt.plot(x, func(x), label=r\"$\\tilde{p}(z)$\")\n",
    "plt.hist(samples, normed=True, alpha=0.2)\n",
    "plt.scatter(samples, np.random.normal(scale=.03, size=(100, 1)), s=5, label=\"samples\")\n",
    "plt.legend(fontsize=15)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
